[2]  1  // Copyright (C) 2002, International Business Machines 

 2  // Corporation and others. All Rights Reserved. 

 3  #if defined(_MSC_VER) 

 4  // Turn off compiler warning about long names 

 5  # pragma warning(disable:4786) 

 6  #endif 

[311]  7  

[325]  8  #include "CbcConfig.h" 

[311]  9  

[2]  10  #include <string> 

 11  //#define CBC_DEBUG 1 

 12  //#define CHECK_CUT_COUNTS 

[135]  13  //#define CHECK_NODE 

[687]  14  //#define CBC_CHECK_BASIS 

[267]  15  #define CBC_WEAK_STRONG 

[2]  16  #include <cassert> 

 17  #include <cfloat> 

 18  #define CUTS 

 19  #include "OsiSolverInterface.hpp" 

[640]  20  #include "OsiChooseVariable.hpp" 

[264]  21  #include "OsiAuxInfo.hpp" 

[222]  22  #include "OsiSolverBranch.hpp" 

[2]  23  #include "CoinWarmStartBasis.hpp" 

[35]  24  #include "CoinTime.hpp" 

[2]  25  #include "CbcModel.hpp" 

 26  #include "CbcNode.hpp" 

[259]  27  #include "CbcStatistics.hpp" 

[271]  28  #include "CbcStrategy.hpp" 

[2]  29  #include "CbcBranchActual.hpp" 

[135]  30  #include "CbcBranchDynamic.hpp" 

[2]  31  #include "OsiRowCut.hpp" 

 32  #include "OsiRowCutDebugger.hpp" 

 33  #include "OsiCuts.hpp" 

 34  #include "CbcCountRowCut.hpp" 

[165]  35  #include "CbcFeasibilityBase.hpp" 

[2]  36  #include "CbcMessage.hpp" 

[311]  37  #ifdef COIN_HAS_CLP 

[2]  38  #include "OsiClpSolverInterface.hpp" 

[122]  39  #include "ClpSimplexOther.hpp" 

[277]  40  #endif 

[2]  41  using namespace std; 

 42  #include "CglCutGenerator.hpp" 

 43  // Default Constructor 

 44  CbcNodeInfo::CbcNodeInfo () 

 45  : 

 46  numberPointingToThis_(0), 

 47  parent_(NULL), 

[912]  48  parentBranch_(NULL), 

[2]  49  owner_(NULL), 

 50  numberCuts_(0), 

[146]  51  nodeNumber_(0), 

[2]  52  cuts_(NULL), 

 53  numberRows_(0), 

[838]  54  numberBranchesLeft_(0), 

 55  active_(7) 

[2]  56  { 

 57  #ifdef CHECK_NODE 

 58  printf("CbcNodeInfo %x Constructor\n",this); 

 59  #endif 

 60  } 

[912]  61  

 62  void 

 63  CbcNodeInfo::setParentBasedData() 

 64  { 

 65  if (parent_) { 

 66  numberRows_ = parent_>numberRows_+parent_>numberCuts_; 

 67  //parent_>increment(); 

 68  if (parent_>owner()) { 

 69  const OsiBranchingObject* br = parent_>owner()>branchingObject(); 

[915]  70  assert(br); 

 71  parentBranch_ = br>clone(); 

[912]  72  } 

 73  } 

 74  } 

 75  

 76  #if 0 

[2]  77  // Constructor given parent 

 78  CbcNodeInfo::CbcNodeInfo (CbcNodeInfo * parent) 

 79  : 

 80  numberPointingToThis_(2), 

 81  parent_(parent), 

[912]  82  parentBranch_(NULL), 

[2]  83  owner_(NULL), 

 84  numberCuts_(0), 

[146]  85  nodeNumber_(0), 

[2]  86  cuts_(NULL), 

 87  numberRows_(0), 

[838]  88  numberBranchesLeft_(2), 

 89  active_(7) 

[2]  90  { 

 91  #ifdef CHECK_NODE 

 92  printf("CbcNodeInfo %x Constructor from parent %x\n",this,parent_); 

 93  #endif 

[912]  94  setParentBasedData(); 

[2]  95  } 

[912]  96  #endif 

 97  

[2]  98  // Copy Constructor 

 99  CbcNodeInfo::CbcNodeInfo (const CbcNodeInfo & rhs) 

 100  : 

 101  numberPointingToThis_(rhs.numberPointingToThis_), 

 102  parent_(rhs.parent_), 

[912]  103  parentBranch_(NULL), 

[2]  104  owner_(rhs.owner_), 

 105  numberCuts_(rhs.numberCuts_), 

[146]  106  nodeNumber_(rhs.nodeNumber_), 

[2]  107  cuts_(NULL), 

 108  numberRows_(rhs.numberRows_), 

[838]  109  numberBranchesLeft_(rhs.numberBranchesLeft_), 

 110  active_(rhs.active_) 

[2]  111  { 

 112  #ifdef CHECK_NODE 

 113  printf("CbcNodeInfo %x Copy constructor\n",this); 

 114  #endif 

 115  if (numberCuts_) { 

 116  cuts_ = new CbcCountRowCut * [numberCuts_]; 

[22]  117  int n=0; 

[2]  118  for (int i=0;i<numberCuts_;i++) { 

 119  CbcCountRowCut * thisCut = rhs.cuts_[i]; 

[22]  120  if (thisCut) { 

 121  // I think this is correct  new one should take priority 

 122  thisCut>setInfo(this,n); 

 123  thisCut>increment(numberBranchesLeft_); 

 124  cuts_[n++] = thisCut; 

 125  } 

[2]  126  } 

[22]  127  numberCuts_=n; 

[2]  128  } 

[912]  129  if (rhs.parentBranch_) { 

 130  parentBranch_ = rhs.parentBranch_>clone(); 

 131  } 

[2]  132  } 

 133  // Constructor given parent and owner 

 134  CbcNodeInfo::CbcNodeInfo (CbcNodeInfo * parent, CbcNode * owner) 

 135  : 

 136  numberPointingToThis_(2), 

 137  parent_(parent), 

[912]  138  parentBranch_(NULL), 

[2]  139  owner_(owner), 

 140  numberCuts_(0), 

[146]  141  nodeNumber_(0), 

[2]  142  cuts_(NULL), 

 143  numberRows_(0), 

[838]  144  numberBranchesLeft_(2), 

 145  active_(7) 

[2]  146  { 

 147  #ifdef CHECK_NODE 

 148  printf("CbcNodeInfo %x Constructor from parent %x\n",this,parent_); 

 149  #endif 

[912]  150  setParentBasedData(); 

[2]  151  } 

 152  

 153  /** 

 154  Take care to detach from the owning CbcNode and decrement the reference 

 155  count in the parent. If this is the last nodeInfo object pointing to the 

 156  parent, make a recursive call to delete the parent. 

 157  */ 

 158  CbcNodeInfo::~CbcNodeInfo() 

 159  { 

 160  #ifdef CHECK_NODE 

 161  printf("CbcNodeInfo %x Destructor parent %x\n",this,parent_); 

 162  #endif 

 163  

 164  assert(!numberPointingToThis_); 

[838]  165  // But there may be some left (max nodes?) 

[640]  166  for (int i=0;i<numberCuts_;i++) { 

[838]  167  if (cuts_[i]) { 

[640]  168  #ifndef GLOBAL_CUTS_JUST_POINTERS 

[838]  169  delete cuts_[i]; 

[640]  170  #else 

[838]  171  if (cuts_[i]>globallyValidAsInteger()!=2) 

 172  delete cuts_[i]; 

[640]  173  #endif 

[838]  174  } 

[640]  175  } 

[2]  176  delete [] cuts_; 

 177  if (owner_) 

 178  owner_>nullNodeInfo(); 

 179  if (parent_) { 

 180  int numberLinks = parent_>decrement(); 

 181  if (!numberLinks) delete parent_; 

 182  } 

[912]  183  delete parentBranch_; 

[2]  184  } 

 185  

 186  

 187  //#define ALLCUTS 

 188  void 

 189  CbcNodeInfo::decrementCuts(int change) 

 190  { 

 191  int i; 

 192  // get rid of all remaining if negative 

 193  int changeThis; 

 194  if (change<0) 

 195  changeThis = numberBranchesLeft_; 

 196  else 

 197  changeThis = change; 

 198  // decrement cut counts 

 199  for (i=0;i<numberCuts_;i++) { 

 200  if (cuts_[i]) { 

 201  int number = cuts_[i]>decrement(changeThis); 

 202  if (!number) { 

 203  //printf("info %x del cut %d %x\n",this,i,cuts_[i]); 

[640]  204  #ifndef GLOBAL_CUTS_JUST_POINTERS 

[2]  205  delete cuts_[i]; 

[640]  206  #else 

 207  if (cuts_[i]>globallyValidAsInteger()!=2) 

 208  delete cuts_[i]; 

 209  #endif 

[2]  210  cuts_[i]=NULL; 

 211  } 

 212  } 

 213  } 

 214  } 

 215  void 

[838]  216  CbcNodeInfo::incrementCuts(int change) 

 217  { 

 218  int i; 

 219  assert (change>0); 

 220  // increment cut counts 

 221  for (i=0;i<numberCuts_;i++) { 

 222  if (cuts_[i]) 

 223  cuts_[i]>increment(change); 

 224  } 

 225  } 

 226  void 

[854]  227  CbcNodeInfo::decrementParentCuts(CbcModel * model,int change) 

[2]  228  { 

 229  if (parent_) { 

 230  // get rid of all remaining if negative 

 231  int changeThis; 

 232  if (change<0) 

 233  changeThis = numberBranchesLeft_; 

 234  else 

 235  changeThis = change; 

 236  int i; 

 237  // Get overestimate of space needed for basis 

[854]  238  CoinWarmStartBasis & dummy = model>workingBasis(); 

[2]  239  dummy.setSize(0,numberRows_+numberCuts_); 

 240  buildRowBasis(dummy); 

 241  /* everything is zero (i.e. free) so we can use to see 

 242  if latest basis */ 

 243  CbcNodeInfo * thisInfo = parent_; 

 244  while (thisInfo) 

 245  thisInfo = thisInfo>buildRowBasis(dummy); 

 246  // decrement cut counts 

 247  thisInfo = parent_; 

 248  int numberRows=numberRows_; 

 249  while (thisInfo) { 

 250  for (i=thisInfo>numberCuts_1;i>=0;i) { 

 251  CoinWarmStartBasis::Status status = dummy.getArtifStatus(numberRows); 

 252  #ifdef ALLCUTS 

 253  status = CoinWarmStartBasis::isFree; 

 254  #endif 

 255  if (thisInfo>cuts_[i]) { 

 256  int number=1; 

 257  if (status!=CoinWarmStartBasis::basic) { 

 258  // tight  drop 1 or 2 

 259  if (change<0) 

 260  number = thisInfo>cuts_[i]>decrement(changeThis); 

 261  else 

 262  number = thisInfo>cuts_[i]>decrement(change); 

 263  } 

 264  if (!number) { 

[640]  265  #ifndef GLOBAL_CUTS_JUST_POINTERS 

[2]  266  delete thisInfo>cuts_[i]; 

[640]  267  #else 

 268  if (thisInfo>cuts_[i]>globallyValidAsInteger()!=2) 

 269  delete thisInfo>cuts_[i]; 

 270  #endif 

[2]  271  thisInfo>cuts_[i]=NULL; 

 272  } 

 273  } 

 274  } 

 275  thisInfo = thisInfo>parent_; 

 276  } 

 277  } 

 278  } 

[854]  279  #if 0 

[2]  280  void 

[854]  281  CbcNodeInfo::incrementParentCuts(CbcModel * model, int change) 

[2]  282  { 

 283  if (parent_) { 

 284  int i; 

 285  // Get overestimate of space needed for basis 

[854]  286  CoinWarmStartBasis & dummy = model>workingBasis(); 

[2]  287  dummy.setSize(0,numberRows_+numberCuts_); 

 288  /* everything is zero (i.e. free) so we can use to see 

 289  if latest basis */ 

 290  buildRowBasis(dummy); 

 291  CbcNodeInfo * thisInfo = parent_; 

 292  while (thisInfo) 

 293  thisInfo = thisInfo>buildRowBasis(dummy); 

 294  // increment cut counts 

 295  thisInfo = parent_; 

 296  int numberRows=numberRows_; 

 297  while (thisInfo) { 

 298  for (i=thisInfo>numberCuts_1;i>=0;i) { 

 299  CoinWarmStartBasis::Status status = dummy.getArtifStatus(numberRows); 

 300  #ifdef ALLCUTS 

 301  status = CoinWarmStartBasis::isFree; 

 302  #endif 

 303  if (thisInfo>cuts_[i]&&status!=CoinWarmStartBasis::basic) { 

 304  thisInfo>cuts_[i]>increment(change); 

 305  } 

 306  } 

 307  thisInfo = thisInfo>parent_; 

 308  } 

 309  } 

 310  } 

[854]  311  #endif 

[2]  312  /* 

 313  Append cuts to the cuts_ array in a nodeInfo. The initial reference count 

 314  is set to numberToBranchOn, which will normally be the number of arms 

 315  defined for the CbcBranchingObject attached to the CbcNode that owns this 

 316  CbcNodeInfo. 

 317  */ 

 318  void 

 319  CbcNodeInfo::addCuts (OsiCuts & cuts, int numberToBranchOn, 

 320  int * whichGenerator) 

 321  { 

 322  int numberCuts = cuts.sizeRowCuts(); 

 323  if (numberCuts) { 

 324  int i; 

 325  if (!numberCuts_) { 

 326  cuts_ = new CbcCountRowCut * [numberCuts]; 

 327  } else { 

 328  CbcCountRowCut ** temp = new CbcCountRowCut * [numberCuts+numberCuts_]; 

 329  memcpy(temp,cuts_,numberCuts_*sizeof(CbcCountRowCut *)); 

 330  delete [] cuts_; 

 331  cuts_ = temp; 

 332  } 

 333  for (i=0;i<numberCuts;i++) { 

 334  CbcCountRowCut * thisCut = new CbcCountRowCut(*cuts.rowCutPtr(i), 

 335  this,numberCuts_); 

 336  thisCut>increment(numberToBranchOn); 

 337  cuts_[numberCuts_++] = thisCut; 

 338  #ifdef CBC_DEBUG 

[300]  339  #if CBC_DEBUG>1 

[2]  340  int n=thisCut>row().getNumElements(); 

 341  printf("Cut %d has %d entries, rhs %g %g =>",i,n,thisCut>lb(), 

 342  thisCut>ub()); 

 343  int j; 

 344  const int * index = thisCut>row().getIndices(); 

 345  const double * element = thisCut>row().getElements(); 

 346  for (j=0;j<n;j++) { 

 347  printf(" (%d,%g)",index[j],element[j]); 

 348  assert(fabs(element[j])>1.00e12); 

 349  } 

 350  printf("\n"); 

[300]  351  #else 

 352  int n=thisCut>row().getNumElements(); 

 353  int j; 

 354  const double * element = thisCut>row().getElements(); 

 355  for (j=0;j<n;j++) { 

 356  assert(fabs(element[j])>1.00e12); 

 357  } 

[2]  358  #endif 

[300]  359  #endif 

[2]  360  } 

 361  } 

 362  } 

 363  

 364  void 

 365  CbcNodeInfo::addCuts(int numberCuts, CbcCountRowCut ** cut, 

 366  int numberToBranchOn) 

 367  { 

 368  if (numberCuts) { 

 369  int i; 

 370  if (!numberCuts_) { 

 371  cuts_ = new CbcCountRowCut * [numberCuts]; 

 372  } else { 

 373  CbcCountRowCut ** temp = new CbcCountRowCut * [numberCuts+numberCuts_]; 

 374  memcpy(temp,cuts_,numberCuts_*sizeof(CbcCountRowCut *)); 

 375  delete [] cuts_; 

 376  cuts_ = temp; 

 377  } 

 378  for (i=0;i<numberCuts;i++) { 

 379  CbcCountRowCut * thisCut = cut[i]; 

 380  thisCut>setInfo(this,numberCuts_); 

 381  //printf("info %x cut %d %x\n",this,i,thisCut); 

 382  thisCut>increment(numberToBranchOn); 

 383  cuts_[numberCuts_++] = thisCut; 

 384  #ifdef CBC_DEBUG 

 385  int n=thisCut>row().getNumElements(); 

 386  #if CBC_DEBUG>1 

 387  printf("Cut %d has %d entries, rhs %g %g =>",i,n,thisCut>lb(), 

 388  thisCut>ub()); 

 389  #endif 

 390  int j; 

 391  #if CBC_DEBUG>1 

 392  const int * index = thisCut>row().getIndices(); 

 393  #endif 

 394  const double * element = thisCut>row().getElements(); 

 395  for (j=0;j<n;j++) { 

 396  #if CBC_DEBUG>1 

 397  printf(" (%d,%g)",index[j],element[j]); 

 398  #endif 

 399  assert(fabs(element[j])>1.00e12); 

 400  } 

 401  printf("\n"); 

 402  #endif 

 403  } 

 404  } 

 405  } 

 406  

 407  // delete cuts 

 408  void 

 409  CbcNodeInfo::deleteCuts(int numberToDelete, CbcCountRowCut ** cuts) 

 410  { 

 411  int i; 

 412  int j; 

 413  int last=1; 

 414  for (i=0;i<numberToDelete;i++) { 

 415  CbcCountRowCut * next = cuts[i]; 

 416  for (j=last+1;j<numberCuts_;j++) { 

 417  if (next==cuts_[j]) 

 418  break; 

 419  } 

 420  if (j==numberCuts_) { 

 421  // start from beginning 

 422  for (j=0;j<last;j++) { 

 423  if (next==cuts_[j]) 

 424  break; 

 425  } 

 426  assert(j<last); 

 427  } 

 428  last=j; 

 429  int number = cuts_[j]>decrement(); 

[640]  430  if (!number) { 

 431  #ifndef GLOBAL_CUTS_JUST_POINTERS 

[2]  432  delete cuts_[j]; 

[640]  433  #else 

 434  if (cuts_[j]>globallyValidAsInteger()!=2) 

 435  delete cuts_[j]; 

 436  #endif 

 437  } 

[2]  438  cuts_[j]=NULL; 

 439  } 

 440  j=0; 

 441  for (i=0;i<numberCuts_;i++) { 

 442  if (cuts_[i]) 

 443  cuts_[j++]=cuts_[i]; 

 444  } 

 445  numberCuts_ = j; 

 446  } 

 447  

 448  // delete cuts 

 449  void 

 450  CbcNodeInfo::deleteCuts(int numberToDelete, int * which) 

 451  { 

 452  int i; 

 453  for (i=0;i<numberToDelete;i++) { 

 454  int iCut=which[i]; 

 455  int number = cuts_[iCut]>decrement(); 

[640]  456  if (!number) { 

 457  #ifndef GLOBAL_CUTS_JUST_POINTERS 

[2]  458  delete cuts_[iCut]; 

[640]  459  #else 

 460  if (cuts_[iCut]>globallyValidAsInteger()!=2) 

 461  delete cuts_[iCut]; 

 462  #endif 

 463  } 

[2]  464  cuts_[iCut]=NULL; 

 465  } 

 466  int j=0; 

 467  for (i=0;i<numberCuts_;i++) { 

 468  if (cuts_[i]) 

 469  cuts_[j++]=cuts_[i]; 

 470  } 

 471  numberCuts_ = j; 

 472  } 

 473  

 474  // Really delete a cut 

 475  void 

 476  CbcNodeInfo::deleteCut(int whichOne) 

 477  { 

 478  assert(whichOne<numberCuts_); 

 479  cuts_[whichOne]=NULL; 

 480  } 

[838]  481  /* Deactivate node information. 

 482  1  bounds 

 483  2  cuts 

 484  4  basis! 

 485  */ 

 486  void 

 487  CbcNodeInfo::deactivate(int mode) 

 488  { 

 489  active_ &= (~mode); 

 490  } 

[2]  491  

 492  CbcFullNodeInfo::CbcFullNodeInfo() : 

 493  CbcNodeInfo(), 

 494  basis_(), 

 495  numberIntegers_(0), 

 496  lower_(NULL), 

 497  upper_(NULL) 

 498  { 

 499  } 

 500  CbcFullNodeInfo::CbcFullNodeInfo(CbcModel * model, 

 501  int numberRowsAtContinuous) : 

[912]  502  CbcNodeInfo(NULL, model>currentNode()) 

[2]  503  { 

 504  OsiSolverInterface * solver = model>solver(); 

 505  numberRows_ = numberRowsAtContinuous; 

 506  numberIntegers_ = model>numberIntegers(); 

 507  int numberColumns = model>getNumCols(); 

 508  lower_ = new double [numberColumns]; 

 509  upper_ = new double [numberColumns]; 

 510  const double * lower = solver>getColLower(); 

 511  const double * upper = solver>getColUpper(); 

 512  int i; 

 513  

 514  for (i=0;i<numberColumns;i++) { 

 515  lower_[i]=lower[i]; 

 516  upper_[i]=upper[i]; 

 517  } 

 518  

 519  basis_ = dynamic_cast<CoinWarmStartBasis*>(solver>getWarmStart()); 

 520  } 

 521  

 522  CbcFullNodeInfo::CbcFullNodeInfo(const CbcFullNodeInfo & rhs) : 

 523  CbcNodeInfo(rhs) 

 524  { 

 525  basis_= dynamic_cast<CoinWarmStartBasis *>(rhs.basis_>clone()) ; 

 526  numberIntegers_=rhs.numberIntegers_; 

 527  lower_=NULL; 

 528  upper_=NULL; 

 529  if (rhs.lower_!=NULL) { 

 530  int numberColumns = basis_>getNumStructural(); 

 531  lower_ = new double [numberColumns]; 

 532  upper_ = new double [numberColumns]; 

 533  assert (upper_!=NULL); 

 534  memcpy(lower_,rhs.lower_,numberColumns*sizeof(double)); 

 535  memcpy(upper_,rhs.upper_,numberColumns*sizeof(double)); 

 536  } 

 537  } 

 538  

 539  CbcNodeInfo * 

 540  CbcFullNodeInfo::clone() const 

 541  { 

 542  return (new CbcFullNodeInfo(*this)); 

 543  } 

 544  

 545  CbcFullNodeInfo::~CbcFullNodeInfo () 

 546  { 

 547  delete basis_ ; 

 548  delete [] lower_; 

 549  delete [] upper_; 

 550  } 

 551  

 552  /* 

 553  The basis supplied as a parameter is deleted and replaced with a new basis 

 554  appropriate for the node, and lower and upper bounds on variables are 

 555  reset according to the stored bounds arrays. Any cuts associated with this 

 556  node are added to the list in addCuts, but not actually added to the 

 557  constraint system in the model. 

 558  

 559  Why pass in a basis at all? The short answer is ``We need the parameter to 

 560  pass out a basis, so might as well use it to pass in the size.'' 

 561  

 562  A longer answer is that in practice we take a memory allocation hit up in 

 563  addCuts1 (the only place applyToModel is called) when we setSize() the 

 564  basis that's passed in. It's immediately tossed here in favour of a clone 

 565  of the basis attached to this nodeInfo. This can probably be fixed, given 

 566  a bit of thought. 

 567  */ 

 568  

 569  void CbcFullNodeInfo::applyToModel (CbcModel *model, 

 570  CoinWarmStartBasis *&basis, 

 571  CbcCountRowCut **addCuts, 

 572  int ¤tNumberCuts) const 

 573  

 574  { OsiSolverInterface *solver = model>solver() ; 

 575  

 576  // branch  do bounds 

[838]  577  assert (active_==7active_==15); 

[2]  578  int i; 

[197]  579  solver>setColLower(lower_); 

 580  solver>setColUpper(upper_); 

[2]  581  int numberColumns = model>getNumCols(); 

 582  // move basis  but make sure size stays 

[264]  583  // for bonmin  should not be needed int numberRows = model>getNumRows(); 

 584  int numberRows=basis>getNumArtificial(); 

[2]  585  delete basis ; 

[264]  586  if (basis_) { 

 587  basis = dynamic_cast<CoinWarmStartBasis *>(basis_>clone()) ; 

 588  basis>resize(numberRows,numberColumns); 

[931]  589  #ifdef CBC_CHECK_BASIS 

 590  std::cout << "Basis (after applying root " <<this<<") "<< std::endl ; 

 591  basis>print() ; 

 592  #endif 

[264]  593  } else { 

 594  // We have a solver without a basis 

 595  basis=NULL; 

 596  } 

[2]  597  for (i=0;i<numberCuts_;i++) 

 598  addCuts[currentNumberCuts+i]= cuts_[i]; 

 599  currentNumberCuts += numberCuts_; 

 600  assert(!parent_); 

 601  return ; 

 602  } 

[838]  603  // Just apply bounds to one variable (1=>infeasible) 

 604  int 

 605  CbcFullNodeInfo::applyBounds(int iColumn, double & lower, double & upper,int force) 

 606  { 

 607  if ((force&&1)==0) { 

 608  if (lower>lower_[iColumn]) 

 609  printf("%d odd lower going from %g to %g\n",iColumn,lower,lower_[iColumn]); 

 610  lower = lower_[iColumn]; 

 611  } else { 

 612  lower_[iColumn]=lower; 

 613  } 

 614  if ((force&&2)==0) { 

 615  if (upper<upper_[iColumn]) 

 616  printf("%d odd upper going from %g to %g\n",iColumn,upper,upper_[iColumn]); 

 617  upper = upper_[iColumn]; 

 618  } else { 

 619  upper_[iColumn]=upper; 

 620  } 

 621  return (upper_[iColumn]>=lower_[iColumn]) ? 0 : 1; 

 622  } 

[2]  623  

 624  /* Builds up row basis backwards (until original model). 

 625  Returns NULL or previous one to apply . 

 626  Depends on Free being 0 and impossible for cuts 

 627  */ 

 628  CbcNodeInfo * 

 629  CbcFullNodeInfo::buildRowBasis(CoinWarmStartBasis & basis ) const 

 630  { 

 631  const unsigned int * saved = 

 632  (const unsigned int *) basis_>getArtificialStatus(); 

 633  unsigned int * now = 

 634  (unsigned int *) basis.getArtificialStatus(); 

 635  int number=basis_>getNumArtificial()>>4;; 

 636  int i; 

 637  for (i=0;i<number;i++) { 

 638  if (!now[i]) 

 639  now[i] = saved[i]; 

 640  } 

 641  return NULL; 

 642  } 

 643  

[819]  644  

[2]  645  // Default constructor 

 646  CbcPartialNodeInfo::CbcPartialNodeInfo() 

 647  

 648  : CbcNodeInfo(), 

 649  basisDiff_(NULL), 

 650  variables_(NULL), 

 651  newBounds_(NULL), 

 652  numberChangedBounds_(0) 

 653  

 654  { /* this space intentionally left blank */ } 

 655  

 656  // Constructor from current state 

 657  CbcPartialNodeInfo::CbcPartialNodeInfo (CbcNodeInfo *parent, CbcNode *owner, 

 658  int numberChangedBounds, 

 659  const int *variables, 

 660  const double *boundChanges, 

 661  const CoinWarmStartDiff *basisDiff) 

 662  : CbcNodeInfo(parent,owner) 

 663  { 

 664  basisDiff_ = basisDiff>clone() ; 

 665  

 666  numberChangedBounds_ = numberChangedBounds; 

[819]  667  int size = numberChangedBounds_*(sizeof(double)+sizeof(int)); 

 668  char * temp = new char [size]; 

 669  newBounds_ = (double *) temp; 

 670  variables_ = (int *) (newBounds_+numberChangedBounds_); 

[2]  671  

 672  int i ; 

 673  for (i=0;i<numberChangedBounds_;i++) { 

 674  variables_[i]=variables[i]; 

 675  newBounds_[i]=boundChanges[i]; 

 676  } 

 677  } 

 678  

 679  CbcPartialNodeInfo::CbcPartialNodeInfo (const CbcPartialNodeInfo & rhs) 

 680  

[912]  681  : CbcNodeInfo(rhs) 

[2]  682  

 683  { basisDiff_ = rhs.basisDiff_>clone() ; 

 684  

 685  numberChangedBounds_ = rhs.numberChangedBounds_; 

[819]  686  int size = numberChangedBounds_*(sizeof(double)+sizeof(int)); 

 687  char * temp = new char [size]; 

 688  newBounds_ = (double *) temp; 

 689  variables_ = (int *) (newBounds_+numberChangedBounds_); 

[2]  690  

 691  int i ; 

 692  for (i=0;i<numberChangedBounds_;i++) { 

 693  variables_[i]=rhs.variables_[i]; 

 694  newBounds_[i]=rhs.newBounds_[i]; 

 695  } 

 696  } 

 697  

 698  CbcNodeInfo * 

 699  CbcPartialNodeInfo::clone() const 

 700  { 

 701  return (new CbcPartialNodeInfo(*this)); 

 702  } 

 703  

 704  

 705  CbcPartialNodeInfo::~CbcPartialNodeInfo () 

 706  { 

 707  delete basisDiff_ ; 

 708  delete [] newBounds_; 

 709  } 

 710  

 711  

 712  /** 

 713  The basis supplied as a parameter is incrementally modified, and lower and 

 714  upper bounds on variables in the model are incrementally modified. Any 

 715  cuts associated with this node are added to the list in addCuts. 

 716  */ 

 717  

 718  void CbcPartialNodeInfo::applyToModel (CbcModel *model, 

 719  CoinWarmStartBasis *&basis, 

 720  CbcCountRowCut **addCuts, 

 721  int ¤tNumberCuts) const 

 722  

 723  { OsiSolverInterface *solver = model>solver(); 

[838]  724  if ((active_&4)!=0) { 

 725  basis>applyDiff(basisDiff_) ; 

[931]  726  #ifdef CBC_CHECK_BASIS 

 727  std::cout << "Basis (after applying " <<this<<") "<< std::endl ; 

 728  basis>print() ; 

 729  #endif 

[838]  730  } 

[2]  731  

 732  // branch  do bounds 

 733  int i; 

[838]  734  if ((active_&1)!=0) { 

 735  for (i=0;i<numberChangedBounds_;i++) { 

 736  int variable = variables_[i]; 

 737  int k = variable&0x3fffffff; 

 738  if ((variable&0x80000000)==0) { 

 739  // lower bound changing 

 740  //#define CBC_PRINT2 

 741  #ifdef CBC_PRINT2 

 742  if(solver>getColLower()[k]!=newBounds_[i]) 

 743  printf("lower change for column %d  from %g to %g\n",k,solver>getColLower()[k],newBounds_[i]); 

 744  #endif 

[640]  745  #ifndef NDEBUG 

[838]  746  if ((variable&0x40000000)==0&&false) { 

 747  double oldValue = solver>getColLower()[k]; 

 748  assert (newBounds_[i]>oldValue1.0e8); 

 749  if (newBounds_[i]<oldValue+1.0e8) 

 750  printf("bad null lower change for column %d  bound %g\n",k,oldValue); 

 751  } 

[640]  752  #endif 

[838]  753  solver>setColLower(k,newBounds_[i]); 

 754  } else { 

 755  // upper bound changing 

 756  #ifdef CBC_PRINT2 

 757  if(solver>getColUpper()[k]!=newBounds_[i]) 

 758  printf("upper change for column %d  from %g to %g\n",k,solver>getColUpper()[k],newBounds_[i]); 

 759  #endif 

[640]  760  #ifndef NDEBUG 

[838]  761  if ((variable&0x40000000)==0&&false) { 

 762  double oldValue = solver>getColUpper()[k]; 

 763  assert (newBounds_[i]<oldValue+1.0e8); 

 764  if (newBounds_[i]>oldValue1.0e8) 

 765  printf("bad null upper change for column %d  bound %g\n",k,oldValue); 

 766  } 

[640]  767  #endif 

[838]  768  solver>setColUpper(k,newBounds_[i]); 

 769  } 

[2]  770  } 

 771  } 

[838]  772  if ((active_&2)!=0) { 

 773  for (i=0;i<numberCuts_;i++) { 

 774  addCuts[currentNumberCuts+i]= cuts_[i]; 

 775  if (cuts_[i]&&model>messageHandler()>logLevel()>4) { 

 776  cuts_[i]>print(); 

 777  } 

[216]  778  } 

[838]  779  

 780  currentNumberCuts += numberCuts_; 

[216]  781  } 

[2]  782  return ; 

 783  } 

[838]  784  // Just apply bounds to one variable (1=>infeasible) 

 785  int 

 786  CbcPartialNodeInfo::applyBounds(int iColumn, double & lower, double & upper,int force) 

 787  { 

 788  // branch  do bounds 

 789  int i; 

 790  int found=0; 

 791  double newLower = COIN_DBL_MAX; 

 792  double newUpper = COIN_DBL_MAX; 

 793  for (i=0;i<numberChangedBounds_;i++) { 

 794  int variable = variables_[i]; 

 795  int k = variable&0x3fffffff; 

 796  if (k==iColumn) { 

 797  if ((variable&0x80000000)==0) { 

 798  // lower bound changing 

 799  found = 1; 

 800  newLower = CoinMax(newLower,newBounds_[i]); 

 801  if ((force&1)==0) { 

 802  if (lower>newBounds_[i]) 

 803  printf("%d odd lower going from %g to %g\n",iColumn,lower,newBounds_[i]); 

 804  lower = newBounds_[i]; 

 805  } else { 

 806  newBounds_[i]=lower; 

 807  variables_[i] = 0x40000000; // say can go odd way 

 808  } 

 809  } else { 

 810  // upper bound changing 

 811  found = 2; 

 812  newUpper = CoinMin(newUpper,newBounds_[i]); 

 813  if ((force&2)==0) { 

 814  if (upper<newBounds_[i]) 

 815  printf("%d odd upper going from %g to %g\n",iColumn,upper,newBounds_[i]); 

 816  upper = newBounds_[i]; 

 817  } else { 

 818  newBounds_[i]=upper; 

 819  variables_[i] = 0x40000000; // say can go odd way 

 820  } 

 821  } 

 822  } 

 823  } 

 824  newLower = CoinMax(newLower,lower); 

 825  newUpper = CoinMin(newUpper,upper); 

 826  int nAdd=0; 

 827  if ((force&2)!=0&&(found&2)==0) { 

 828  // need to add new upper 

 829  nAdd++; 

 830  } 

 831  if ((force&1)!=0&&(found&1)==0) { 

 832  // need to add new lower 

 833  nAdd++; 

 834  } 

 835  if (nAdd) { 

 836  int size = (numberChangedBounds_+nAdd)*(sizeof(double)+sizeof(int)); 

 837  char * temp = new char [size]; 

 838  double * newBounds = (double *) temp; 

 839  int * variables = (int *) (newBounds+numberChangedBounds_+nAdd); 

[2]  840  

[838]  841  int i ; 

 842  for (i=0;i<numberChangedBounds_;i++) { 

 843  variables[i]=variables_[i]; 

 844  newBounds[i]=newBounds_[i]; 

 845  } 

 846  delete [] newBounds_; 

 847  newBounds_ = newBounds; 

 848  variables_ = variables; 

 849  if ((force&2)!=0&&(found&2)==0) { 

 850  // need to add new upper 

 851  int variable = iColumn  0x80000000; 

 852  variables_[numberChangedBounds_]=variable; 

 853  newBounds_[numberChangedBounds_++]=newUpper; 

 854  } 

 855  if ((force&1)!=0&&(found&1)==0) { 

 856  // need to add new lower 

 857  int variable = iColumn; 

 858  variables_[numberChangedBounds_]=variable; 

 859  newBounds_[numberChangedBounds_++]=newLower; 

 860  } 

 861  } 

 862  

 863  return (newUpper>=newLower) ? 0 : 1; 

 864  } 

 865  

[2]  866  /* Builds up row basis backwards (until original model). 

 867  Returns NULL or previous one to apply . 

 868  Depends on Free being 0 and impossible for cuts 

 869  */ 

 870  

 871  CbcNodeInfo * 

 872  CbcPartialNodeInfo::buildRowBasis(CoinWarmStartBasis & basis ) const 

 873  

 874  { basis.applyDiff(basisDiff_) ; 

 875  

 876  return parent_ ; } 

 877  

 878  CbcNode::CbcNode() : 

 879  nodeInfo_(NULL), 

 880  objectiveValue_(1.0e100), 

 881  guessedObjectiveValue_(1.0e100), 

[460]  882  sumInfeasibilities_(0.0), 

[2]  883  branch_(NULL), 

 884  depth_(1), 

[838]  885  numberUnsatisfied_(0), 

 886  nodeNumber_(1), 

 887  state_(0) 

[2]  888  { 

 889  #ifdef CHECK_NODE 

 890  printf("CbcNode %x Constructor\n",this); 

 891  #endif 

 892  } 

[838]  893  // Print 

 894  void 

 895  CbcNode::print() const 

 896  { 

 897  printf("number %d obj %g depth %d sumun %g nunsat %d state %d\n", 

 898  nodeNumber_,objectiveValue_,depth_,sumInfeasibilities_,numberUnsatisfied_,state_); 

 899  } 

[2]  900  CbcNode::CbcNode(CbcModel * model, 

 901  CbcNode * lastNode) : 

 902  nodeInfo_(NULL), 

 903  objectiveValue_(1.0e100), 

 904  guessedObjectiveValue_(1.0e100), 

[460]  905  sumInfeasibilities_(0.0), 

[2]  906  branch_(NULL), 

 907  depth_(1), 

[838]  908  numberUnsatisfied_(0), 

 909  nodeNumber_(1), 

 910  state_(0) 

[2]  911  { 

 912  #ifdef CHECK_NODE 

 913  printf("CbcNode %x Constructor from model\n",this); 

 914  #endif 

[268]  915  model>setObjectiveValue(this,lastNode); 

[2]  916  

[838]  917  if (lastNode) { 

 918  if (lastNode>nodeInfo_) { 

 919  lastNode>nodeInfo_>increment(); 

 920  } 

 921  } 

 922  nodeNumber_= model>getNodeCount(); 

[2]  923  } 

 924  

[640]  925  #define CBC_NEW_CREATEINFO 

 926  #ifdef CBC_NEW_CREATEINFO 

[2]  927  

[640]  928  /* 

 929  New createInfo, with basis manipulation hidden inside mergeBasis. Allows 

 930  solvers to override and carry over all information from one basis to 

 931  another. 

 932  */ 

 933  

[2]  934  void 

 935  CbcNode::createInfo (CbcModel *model, 

 936  CbcNode *lastNode, 

 937  const CoinWarmStartBasis *lastws, 

 938  const double *lastLower, const double *lastUpper, 

[640]  939  int numberOldActiveCuts, int numberNewCuts) 

 940  

 941  { OsiSolverInterface *solver = model>solver(); 

 942  CbcStrategy *strategy = model>strategy(); 

 943  /* 

 944  The root  no parent. Create full basis and bounds information. 

 945  */ 

 946  if (!lastNode) 

 947  { 

 948  if (!strategy) 

 949  nodeInfo_=new CbcFullNodeInfo(model,solver>getNumRows()); 

 950  else 

 951  nodeInfo_ = strategy>fullNodeInfo(model,solver>getNumRows()); 

 952  } else { 

 953  /* 

 954  Not the root. Create an edit from the parent's basis & bound information. 

 955  This is not quite as straightforward as it seems. We need to reintroduce 

 956  cuts we may have dropped out of the basis, in the correct position, because 

 957  this whole process is strictly positional. Start by grabbing the current 

 958  basis. 

 959  */ 

[854]  960  bool mustDeleteBasis; 

[640]  961  const CoinWarmStartBasis *ws = 

[854]  962  dynamic_cast<const CoinWarmStartBasis*>(solver>getPointerToWarmStart(mustDeleteBasis)); 

[640]  963  assert(ws!=NULL); // make sure not volume 

 964  //int numberArtificials = lastws>getNumArtificial(); 

 965  int numberColumns = solver>getNumCols(); 

 966  int numberRowsAtContinuous = model>numberRowsAtContinuous(); 

 967  int currentNumberCuts = model>currentNumberCuts(); 

 968  # ifdef CBC_CHECK_BASIS 

 969  std::cout 

 970  << "Before expansion: orig " << numberRowsAtContinuous 

 971  << ", old " << numberOldActiveCuts 

 972  << ", new " << numberNewCuts 

 973  << ", current " << currentNumberCuts << "." << std::endl ; 

 974  ws>print(); 

 975  # endif 

 976  /* 

 977  Clone the basis and resize it to hold the structural constraints, plus 

 978  all the cuts: old cuts, both active and inactive (currentNumberCuts), 

 979  and new cuts (numberNewCuts). This will become the expanded basis. 

 980  */ 

 981  CoinWarmStartBasis *expanded = 

 982  dynamic_cast<CoinWarmStartBasis *>(ws>clone()) ; 

 983  int iCompact = numberRowsAtContinuous+numberOldActiveCuts+numberNewCuts ; 

 984  // int nPartial = numberRowsAtContinuous+currentNumberCuts; 

 985  int iFull = numberRowsAtContinuous+currentNumberCuts+numberNewCuts; 

 986  // int maxBasisLength = ((iFull+15)>>4)+((numberColumns+15)>>4); 

 987  // printf("l %d full %d\n",maxBasisLength,iFull); 

 988  expanded>resize(iFull,numberColumns); 

 989  # ifdef CBC_CHECK_BASIS 

 990  std::cout 

 991  << "\tFull basis " << iFull << " rows, " 

 992  << numberColumns << " columns; compact " 

 993  << iCompact << " rows." << std::endl ; 

 994  # endif 

 995  /* 

 996  Now flesh out the expanded basis. The clone already has the 

 997  correct status information for the variables and for the structural 

 998  (numberRowsAtContinuous) constraints. Any indices beyond nPartial must be 

 999  cuts created while processing this node  they can be copied en bloc 

 1000  into the correct position in the expanded basis. The space reserved for 

 1001  xferRows is a gross overestimate. 

 1002  */ 

 1003  CoinWarmStartBasis::XferVec xferRows ; 

 1004  xferRows.reserve(iFullnumberRowsAtContinuous+1) ; 

 1005  if (numberNewCuts) { 

 1006  xferRows.push_back( 

 1007  CoinWarmStartBasis::XferEntry(iCompactnumberNewCuts, 

 1008  iFullnumberNewCuts,numberNewCuts)) ; 

 1009  } 

 1010  /* 

 1011  From nPartial down, record the entries we want to copy from the current 

 1012  basis (the entries for the active cuts; nonzero in the list returned 

 1013  by addedCuts). Fill the expanded basis with entries showing a status of 

 1014  basic for the deactivated (loose) cuts. 

 1015  */ 

 1016  CbcCountRowCut **cut = model>addedCuts(); 

 1017  iFull = (numberNewCuts+1) ; 

 1018  iCompact = (numberNewCuts+1) ; 

 1019  int runLen = 0 ; 

 1020  CoinWarmStartBasis::XferEntry entry(1,1,1) ; 

 1021  while (iFull >= numberRowsAtContinuous) { 

 1022  for ( ; iFull >= numberRowsAtContinuous && 

 1023  cut[iFullnumberRowsAtContinuous] ; iFull) 

 1024  runLen++ ; 

 1025  if (runLen) { 

 1026  iCompact = runLen ; 

 1027  entry.first = iCompact+1 ; 

 1028  entry.second = iFull+1 ; 

 1029  entry.third = runLen ; 

 1030  runLen = 0 ; 

 1031  xferRows.push_back(entry) ; 

 1032  } 

 1033  for ( ; iFull >= numberRowsAtContinuous && 

 1034  !cut[iFullnumberRowsAtContinuous] ; iFull) 

 1035  expanded>setArtifStatus(iFull,CoinWarmStartBasis::basic); 

 1036  } 

 1037  /* 

 1038  Finally, call mergeBasis to copy over entries from the current basis to 

 1039  the expanded basis. Since we cloned the expanded basis from the active basis 

 1040  and haven't changed the number of variables, only row status entries need 

 1041  to be copied. 

 1042  */ 

 1043  expanded>mergeBasis(ws,&xferRows,0) ; 

 1044  

 1045  #ifdef CBC_CHECK_BASIS 

 1046  std::cout << "Expanded basis:" << std::endl ; 

 1047  expanded>print() ; 

 1048  std::cout << "Diffing against:" << std::endl ; 

 1049  lastws>print() ; 

 1050  #endif 

[931]  1051  assert (expanded>getNumArtificial()>=lastws>getNumArtificial()); 

[1060]  1052  #ifdef CLP_INVESTIGATE 

[931]  1053  if (!expanded>fullBasis()) { 

 1054  int iFull = numberRowsAtContinuous+currentNumberCuts+numberNewCuts; 

 1055  printf("cont %d old %d new %d current %d full inc %d full %d\n", 

 1056  numberRowsAtContinuous,numberOldActiveCuts,numberNewCuts, 

 1057  currentNumberCuts,iFull,iFullnumberNewCuts); 

 1058  } 

[1060]  1059  #endif 

[931]  1060  assert (lastws>fullBasis()); 

[640]  1061  

 1062  /* 

 1063  Now that we have two bases in proper positional correspondence, creating 

 1064  the actual diff is dead easy. 

 1065  

 1066  Note that we're going to compare the expanded basis here to the stripped 

 1067  basis (lastws) produced by addCuts. It doesn't affect the correctness (the 

 1068  diff process has no knowledge of the meaning of an entry) but it does 

 1069  mean that we'll always generate a whack of diff entries because the expanded 

 1070  basis is considerably larger than the stripped basis. 

 1071  */ 

 1072  CoinWarmStartDiff *basisDiff = expanded>generateDiff(lastws) ; 

 1073  /* 

 1074  Diff the bound vectors. It's assumed the number of structural variables 

 1075  is not changing. For branching objects that change bounds on integer 

 1076  variables, we should see at least one bound change as a consequence 

 1077  of applying the branch that generated this subproblem from its parent. 

 1078  This need not hold for other types of branching objects (hyperplane 

 1079  branches, for example). 

 1080  */ 

 1081  const double * lower = solver>getColLower(); 

 1082  const double * upper = solver>getColUpper(); 

 1083  

 1084  double *boundChanges = new double [2*numberColumns] ; 

 1085  int *variables = new int [2*numberColumns] ; 

 1086  int numberChangedBounds=0; 

 1087  

 1088  int i; 

 1089  for (i=0;i<numberColumns;i++) { 

 1090  if (lower[i]!=lastLower[i]) { 

 1091  variables[numberChangedBounds]=i; 

 1092  boundChanges[numberChangedBounds++]=lower[i]; 

 1093  } 

 1094  if (upper[i]!=lastUpper[i]) { 

 1095  variables[numberChangedBounds]=i0x80000000; 

 1096  boundChanges[numberChangedBounds++]=upper[i]; 

 1097  } 

 1098  #ifdef CBC_DEBUG 

 1099  if (lower[i] != lastLower[i]) { 

 1100  std::cout 

 1101  << "lower on " << i << " changed from " 

 1102  << lastLower[i] << " to " << lower[i] << std::endl ; 

 1103  } 

 1104  if (upper[i] != lastUpper[i]) { 

 1105  std::cout 

 1106  << "upper on " << i << " changed from " 

 1107  << lastUpper[i] << " to " << upper[i] << std::endl ; 

 1108  } 

 1109  #endif 

 1110  } 

 1111  #ifdef CBC_DEBUG 

 1112  std::cout << numberChangedBounds << " changed bounds." << std::endl ; 

 1113  #endif 

 1114  //if (lastNode>branchingObject()>boundBranch()) 

 1115  //assert (numberChangedBounds); 

 1116  /* 

 1117  Hand the lot over to the CbcPartialNodeInfo constructor, then clean up and 

 1118  return. 

 1119  */ 

 1120  if (!strategy) 

 1121  nodeInfo_ = 

 1122  new CbcPartialNodeInfo(lastNode>nodeInfo_,this,numberChangedBounds, 

 1123  variables,boundChanges,basisDiff) ; 

 1124  else 

 1125  nodeInfo_ = 

 1126  strategy>partialNodeInfo(model,lastNode>nodeInfo_,this, 

 1127  numberChangedBounds,variables,boundChanges, 

 1128  basisDiff) ; 

 1129  delete basisDiff ; 

 1130  delete [] boundChanges; 

 1131  delete [] variables; 

 1132  delete expanded ; 

[854]  1133  if (mustDeleteBasis) 

 1134  delete ws; 

[640]  1135  } 

 1136  // Set node number 

 1137  nodeInfo_>setNodeNumber(model>getNodeCount2()); 

[838]  1138  state_ = 2; // say active 

[640]  1139  } 

 1140  

 1141  #else // CBC_NEW_CREATEINFO 

 1142  

 1143  /* 

 1144  Original createInfo, with bare manipulation of basis vectors. Fails if solver 

 1145  maintains additional information in basis. 

 1146  */ 

 1147  

 1148  void 

 1149  CbcNode::createInfo (CbcModel *model, 

 1150  CbcNode *lastNode, 

 1151  const CoinWarmStartBasis *lastws, 

 1152  const double *lastLower, const double *lastUpper, 

[2]  1153  int numberOldActiveCuts,int numberNewCuts) 

 1154  { OsiSolverInterface * solver = model>solver(); 

[271]  1155  CbcStrategy * strategy = model>strategy(); 

[2]  1156  /* 

 1157  The root  no parent. Create full basis and bounds information. 

 1158  */ 

 1159  if (!lastNode) 

[271]  1160  { 

 1161  if (!strategy) 

 1162  nodeInfo_=new CbcFullNodeInfo(model,solver>getNumRows()); 

 1163  else 

 1164  nodeInfo_ = strategy>fullNodeInfo(model,solver>getNumRows()); 

 1165  } 

[2]  1166  /* 

 1167  Not the root. Create an edit from the parent's basis & bound information. 

 1168  This is not quite as straightforward as it seems. We need to reintroduce 

 1169  cuts we may have dropped out of the basis, in the correct position, because 

 1170  this whole process is strictly positional. Start by grabbing the current 

 1171  basis. 

 1172  */ 

 1173  else 

[854]  1174  { 

 1175  bool mustDeleteBasis; 

 1176  const CoinWarmStartBasis* ws = 

 1177  dynamic_cast<const CoinWarmStartBasis*>(solver>getPointerToWarmStart(mustDeleteBasis)); 

[2]  1178  assert(ws!=NULL); // make sure not volume 

 1179  //int numberArtificials = lastws>getNumArtificial(); 

 1180  int numberColumns = solver>getNumCols(); 

 1181  

 1182  const double * lower = solver>getColLower(); 

 1183  const double * upper = solver>getColUpper(); 

 1184  

 1185  int i; 

 1186  /* 

 1187  Create a clone and resize it to hold all the structural constraints, plus 

 1188  all the cuts: old cuts, both active and inactive (currentNumberCuts), and 

 1189  new cuts (numberNewCuts). 

 1190  

 1191  TODO: You'd think that the set of constraints (logicals) in the expanded 

 1192  basis should match the set represented in lastws. At least, that's 

 1193  what I thought. But at the point I first looked hard at this bit of 

 1194  code, it turned out that lastws was the stripped basis produced at 

 1195  the end of addCuts(), rather than the raw basis handed back by 

 1196  addCuts1(). The expanded basis here is equivalent to the raw basis of 

 1197  addCuts1(). I said ``whoa, that's not good, I must have introduced a 

 1198  bug'' and went back to John's code to see where I'd gone wrong. 

 1199  And discovered the same `error' in his code. 

 1200  

 1201  After a bit of thought, my conclusion is that correctness is not 

 1202  affected by whether lastws is the stripped or raw basis. The diffs 

 1203  have no semantics  just a set of changes that need to be made 

 1204  to convert lastws into expanded. I think the only effect is that we 

 1205  store a lot more diffs (everything in expanded that's not covered by 

 1206  the stripped basis). But I need to give this more thought. There 

 1207  may well be some subtle error cases. 

 1208  

 1209  In the mean time, I've twiddled addCuts() to set lastws to the raw 

 1210  basis. Makes me (Lou) less nervous to compare apples to apples. 

 1211  */ 

 1212  CoinWarmStartBasis *expanded = 

 1213  dynamic_cast<CoinWarmStartBasis *>(ws>clone()) ; 

 1214  int numberRowsAtContinuous = model>numberRowsAtContinuous(); 

 1215  int iFull = numberRowsAtContinuous+model>currentNumberCuts()+ 

 1216  numberNewCuts; 

 1217  //int numberArtificialsNow = iFull; 

 1218  //int maxBasisLength = ((iFull+15)>>4)+((numberColumns+15)>>4); 

 1219  //printf("l %d full %d\n",maxBasisLength,iFull); 

[264]  1220  if (expanded) 

 1221  expanded>resize(iFull,numberColumns); 

[300]  1222  #ifdef CBC_CHECK_BASIS 

[2]  1223  printf("Before expansion: orig %d, old %d, new %d, current %d\n", 

 1224  numberRowsAtContinuous,numberOldActiveCuts,numberNewCuts, 

 1225  model>currentNumberCuts()) ; 

 1226  ws>print(); 

 1227  #endif 

 1228  /* 

 1229  Now fill in the expanded basis. Any indices beyond nPartial must 

 1230  be cuts created while processing this node  they can be copied directly 

 1231  into the expanded basis. From nPartial down, pull the status of active cuts 

 1232  from ws, interleaving with a B entry for the deactivated (loose) cuts. 

 1233  */ 

 1234  int numberDropped = model>currentNumberCuts()numberOldActiveCuts; 

 1235  int iCompact=iFullnumberDropped; 

 1236  CbcCountRowCut ** cut = model>addedCuts(); 

 1237  int nPartial = model>currentNumberCuts()+numberRowsAtContinuous; 

 1238  iFull; 

 1239  for (;iFull>=nPartial;iFull) { 

 1240  CoinWarmStartBasis::Status status = ws>getArtifStatus(iCompact); 

[13]  1241  //assert (status != CoinWarmStartBasis::basic); // may be permanent cut 

[2]  1242  expanded>setArtifStatus(iFull,status); 

 1243  } 

 1244  for (;iFull>=numberRowsAtContinuous;iFull) { 

 1245  if (cut[iFullnumberRowsAtContinuous]) { 

 1246  CoinWarmStartBasis::Status status = ws>getArtifStatus(iCompact); 

 1247  // If no cut generator being used then we may have basic variables 

[34]  1248  //if (model>getMaximumCutPasses()&& 

 1249  // status == CoinWarmStartBasis::basic) 

 1250  //printf("cut basic\n"); 

[2]  1251  expanded>setArtifStatus(iFull,status); 

 1252  } else { 

 1253  expanded>setArtifStatus(iFull,CoinWarmStartBasis::basic); 

 1254  } 

 1255  } 

[300]  1256  #ifdef CBC_CHECK_BASIS 

[2]  1257  printf("Expanded basis\n"); 

 1258  expanded>print() ; 

 1259  printf("Diffing against\n") ; 

 1260  lastws>print() ; 

 1261  #endif 

 1262  /* 

 1263  Now that we have two bases in proper positional correspondence, creating 

 1264  the actual diff is dead easy. 

 1265  */ 

 1266  

 1267  CoinWarmStartDiff *basisDiff = expanded>generateDiff(lastws) ; 

 1268  /* 

 1269  Diff the bound vectors. It's assumed the number of structural variables is 

 1270  not changing. Assuming that branching objects all involve integer variables, 

 1271  we should see at least one bound change as a consequence of processing this 

 1272  subproblem. Different types of branching objects could break this assertion. 

[150]  1273  Not true at all  we have not applied current branch  JJF. 

[2]  1274  */ 

 1275  double *boundChanges = new double [2*numberColumns] ; 

 1276  int *variables = new int [2*numberColumns] ; 

 1277  int numberChangedBounds=0; 

 1278  for (i=0;i<numberColumns;i++) { 

 1279  if (lower[i]!=lastLower[i]) { 

 1280  variables[numberChangedBounds]=i; 

 1281  boundChanges[numberChangedBounds++]=lower[i]; 

 1282  } 

 1283  if (upper[i]!=lastUpper[i]) { 

 1284  variables[numberChangedBounds]=i0x80000000; 

 1285  boundChanges[numberChangedBounds++]=upper[i]; 

 1286  } 

 1287  #ifdef CBC_DEBUG 

 1288  if (lower[i]!=lastLower[i]) 

 1289  printf("lower on %d changed from %g to %g\n", 

 1290  i,lastLower[i],lower[i]); 

 1291  if (upper[i]!=lastUpper[i]) 

 1292  printf("upper on %d changed from %g to %g\n", 

 1293  i,lastUpper[i],upper[i]); 

 1294  #endif 

 1295  } 

 1296  #ifdef CBC_DEBUG 

 1297  printf("%d changed bounds\n",numberChangedBounds) ; 

 1298  #endif 

[150]  1299  //if (lastNode>branchingObject()>boundBranch()) 

 1300  //assert (numberChangedBounds); 

[2]  1301  /* 

 1302  Hand the lot over to the CbcPartialNodeInfo constructor, then clean up and 

 1303  return. 

 1304  */ 

[271]  1305  if (!strategy) 

 1306  nodeInfo_ = 

 1307  new CbcPartialNodeInfo(lastNode>nodeInfo_,this,numberChangedBounds, 

 1308  variables,boundChanges,basisDiff) ; 

 1309  else 

 1310  nodeInfo_ = strategy>partialNodeInfo(model, lastNode>nodeInfo_,this,numberChangedBounds, 

 1311  variables,boundChanges,basisDiff) ; 

[2]  1312  delete basisDiff ; 

 1313  delete [] boundChanges; 

 1314  delete [] variables; 

 1315  delete expanded ; 

[854]  1316  if (mustDeleteBasis) 

 1317  delete ws; 

[2]  1318  } 

[146]  1319  // Set node number 

[154]  1320  nodeInfo_>setNodeNumber(model>getNodeCount2()); 

[838]  1321  state_ = 2; // say active 

[2]  1322  } 

 1323  

[640]  1324  #endif // CBC_NEW_CREATEINFO 

 1325  

[2]  1326  /* 

 1327  The routine scans through the object list of the model looking for objects 

 1328  that indicate infeasibility. It tests each object using strong branching 

 1329  and selects the one with the least objective degradation. A corresponding 

 1330  branching object is left attached to lastNode. 

 1331  

 1332  If strong branching is disabled, a candidate object is chosen essentially 

 1333  at random (whatever object ends up in pos'n 0 of the candidate array). 

 1334  

 1335  If a branching candidate is found to be monotone, bounds are set to fix the 

 1336  variable and the routine immediately returns (the caller is expected to 

 1337  reoptimize). 

 1338  

 1339  If a branching candidate is found to result in infeasibility in both 

 1340  directions, the routine immediately returns an indication of infeasibility. 

 1341  

 1342  Returns: 0 both branch directions are feasible 

 1343  1 branching variable is monotone 

 1344  2 infeasible 

 1345  

 1346  Original comments: 

 1347  Here could go cuts etc etc 

 1348  For now just fix on objective from strong branching. 

 1349  */ 

 1350  

[48]  1351  int CbcNode::chooseBranch (CbcModel *model, CbcNode *lastNode,int numberPassesLeft) 

[2]  1352  

 1353  { if (lastNode) 

 1354  depth_ = lastNode>depth_+1; 

 1355  else 

 1356  depth_ = 0; 

 1357  delete branch_; 

 1358  branch_=NULL; 

 1359  OsiSolverInterface * solver = model>solver(); 

 1360  double saveObjectiveValue = solver>getObjValue(); 

[268]  1361  double objectiveValue = CoinMax(solver>getObjSense()*saveObjectiveValue,objectiveValue_); 

[2]  1362  const double * lower = solver>getColLower(); 

 1363  const double * upper = solver>getColUpper(); 

[165]  1364  // See what user thinks 

 1365  int anyAction=model>problemFeasibility()>feasible(model,0); 

 1366  if (anyAction) { 

 1367  // will return 2 if infeasible , 0 if treat as integer 

 1368  return anyAction1; 

 1369  } 

[2]  1370  double integerTolerance = 

 1371  model>getDblParam(CbcModel::CbcIntegerTolerance); 

[640]  1372  // point to useful information 

 1373  OsiBranchingInformation usefulInfo = model>usefulInformation(); 

 1374  // and modify 

 1375  usefulInfo.depth_=depth_; 

[2]  1376  int i; 

 1377  bool beforeSolution = model>getSolutionCount()==0; 

 1378  int numberStrong=model>numberStrong(); 

[231]  1379  // switch off strong if hotstart 

 1380  if (model>hotstartSolution()) 

 1381  numberStrong=0; 

[259]  1382  int numberStrongDone=0; 

 1383  int numberUnfinished=0; 

 1384  int numberStrongInfeasible=0; 

 1385  int numberStrongIterations=0; 

[122]  1386  int saveNumberStrong=numberStrong; 

[2]  1387  int numberObjects = model>numberObjects(); 

[104]  1388  bool checkFeasibility = numberObjects>model>numberIntegers(); 

[854]  1389  int maximumStrong = CoinMax(CoinMin(numberStrong,numberObjects),1); 

[2]  1390  int numberColumns = model>getNumCols(); 

 1391  double * saveUpper = new double[numberColumns]; 

 1392  double * saveLower = new double[numberColumns]; 

 1393  

 1394  // Save solution in case heuristics need good solution later 

[122]  1395  

[2]  1396  double * saveSolution = new double[numberColumns]; 

 1397  memcpy(saveSolution,solver>getColSolution(),numberColumns*sizeof(double)); 

[135]  1398  model>reserveCurrentSolution(saveSolution); 

[122]  1399  /* 

 1400  Get a branching decision object. Use the default decision criteria unless 

 1401  the user has loaded a decision method into the model. 

 1402  */ 

[2]  1403  CbcBranchDecision *decision = model>branchingMethod(); 

[640]  1404  CbcDynamicPseudoCostBranchingObject * dynamicBranchingObject = 

 1405  dynamic_cast<CbcDynamicPseudoCostBranchingObject *>(decision); 

 1406  if (!decisiondynamicBranchingObject) 

[2]  1407  decision = new CbcBranchDefaultDecision(); 

[640]  1408  decision>initialize(model); 

[135]  1409  CbcStrongInfo * choice = new CbcStrongInfo[maximumStrong]; 

[2]  1410  for (i=0;i<numberColumns;i++) { 

 1411  saveLower[i] = lower[i]; 

 1412  saveUpper[i] = upper[i]; 

 1413  } 

[122]  1414  // May go round twice if strong branching fixes all local candidates 

 1415  bool finished=false; 

[2]  1416  double estimatedDegradation=0.0; 

[122]  1417  while(!finished) { 

 1418  finished=true; 

[44]  1419  // Some objects may compute an estimate of best solution from here 

 1420  estimatedDegradation=0.0; 

[197]  1421  //int numberIntegerInfeasibilities=0; // without odd ones 

[259]  1422  numberStrongDone=0; 

 1423  numberUnfinished=0; 

 1424  numberStrongInfeasible=0; 

 1425  numberStrongIterations=0; 

[122]  1426  

 1427  // We may go round this loop twice (only if we think we have solution) 

 1428  for (int iPass=0;iPass<2;iPass++) { 

[44]  1429  

[122]  1430  // compute current state 

[197]  1431  //int numberObjectInfeasibilities; // just odd ones 

 1432  //model>feasibleSolution( 

 1433  // numberIntegerInfeasibilities, 

 1434  // numberObjectInfeasibilities); 

[231]  1435  const double * hotstartSolution = model>hotstartSolution(); 

 1436  const int * hotstartPriorities = model>hotstartPriorities(); 

[44]  1437  

[122]  1438  // Some objects may compute an estimate of best solution from here 

 1439  estimatedDegradation=0.0; 

 1440  numberUnsatisfied_ = 0; 

[460]  1441  // initialize sum of "infeasibilities" 

 1442  sumInfeasibilities_ = 0.0; 

[702]  1443  int bestPriority=COIN_INT_MAX; 

[122]  1444  /* 

 1445  Scan for branching objects that indicate infeasibility. Choose the best 

 1446  maximumStrong candidates, using priority as the first criteria, then 

 1447  integer infeasibility. 

 1448  

 1449  The algorithm is to fill the choice array with a set of good candidates (by 

 1450  infeasibility) with priority bestPriority. Finding a candidate with 

 1451  priority better (less) than bestPriority flushes the choice array. (This 

 1452  serves as initialization when the first candidate is found.) 

 1453  

 1454  A new candidate is added to choices only if its infeasibility exceeds the 

 1455  current max infeasibility (mostAway). When a candidate is added, it 

 1456  replaces the candidate with the smallest infeasibility (tracked by 

 1457  iSmallest). 

 1458  */ 

 1459  int iSmallest = 0; 

 1460  double mostAway = 1.0e100; 

[135]  1461  for (i = 0 ; i < maximumStrong ; i++) 

 1462  choice[i].possibleBranch = NULL ; 

[122]  1463  numberStrong=0; 

[307]  1464  bool canDoOneHot=false; 

[122]  1465  for (i=0;i<numberObjects;i++) { 

[640]  1466  OsiObject * object = model>modifiableObject(i); 

[122]  1467  int preferredWay; 

[640]  1468  double infeasibility = object>infeasibility(&usefulInfo,preferredWay); 

[122]  1469  int priorityLevel = object>priority(); 

[231]  1470  if (hotstartSolution) { 

[122]  1471  // we are doing hot start 

 1472  const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object); 

 1473  if (thisOne) { 

[640]  1474  int iColumn = thisOne>columnNumber(); 

[307]  1475  bool canDoThisHot=true; 

[231]  1476  double targetValue = hotstartSolution[iColumn]; 

[122]  1477  if (saveUpper[iColumn]>saveLower[iColumn]) { 

 1478  double value = saveSolution[iColumn]; 

[231]  1479  if (hotstartPriorities) 

 1480  priorityLevel=hotstartPriorities[iColumn]; 

[122]  1481  //double originalLower = thisOne>originalLower(); 

 1482  //double originalUpper = thisOne>originalUpper(); 

 1483  // switch off if not possible 

 1484  if (targetValue>=saveLower[iColumn]&&targetValue<=saveUpper[iColumn]) { 

[231]  1485  /* priority outranks rest always if negative 

[122]  1486  otherwise can be downgraded if at correct level. 

[231]  1487  Infeasibility may be increased to choose 1.0 values first. 

 1488  choose one near wanted value 

[122]  1489  */ 

 1490  if (fabs(valuetargetValue)>integerTolerance) { 

[231]  1491  infeasibility = 1.0fabs(valuetargetValue); 

 1492  if (targetValue==1.0) 

 1493  infeasibility += 1.0; 

[122]  1494  if (value>targetValue) { 

 1495  preferredWay=1; 

 1496  } else { 

 1497  preferredWay=1; 

 1498  } 

[231]  1499  priorityLevel = CoinAbs(priorityLevel); 

 1500  } else if (priorityLevel<0) { 

 1501  priorityLevel = CoinAbs(priorityLevel); 

[122]  1502  if (targetValue==saveLower[iColumn]) { 

[231]  1503  infeasibility = integerTolerance+1.0e12; 

[122]  1504  preferredWay=1; 

 1505  } else if (targetValue==saveUpper[iColumn]) { 

[231]  1506  infeasibility = integerTolerance+1.0e12; 

[122]  1507  preferredWay=1; 

 1508  } else { 

[231]  1509  // can't 

 1510  priorityLevel += 10000000; 

[307]  1511  canDoThisHot=false; 

[122]  1512  } 

[44]  1513  } else { 

[122]  1514  priorityLevel += 10000000; 

[307]  1515  canDoThisHot=false; 

[44]  1516  } 

 1517  } else { 

[122]  1518  // switch off if not possible 

[307]  1519  canDoThisHot=false; 

[44]  1520  } 

[307]  1521  if (canDoThisHot) 

 1522  canDoOneHot=true; 

[231]  1523  } else if (targetValue<saveLower[iColumn]targetValue>saveUpper[iColumn]) { 

[44]  1524  } 

[231]  1525  } else { 

 1526  priorityLevel += 10000000; 

[44]  1527  } 

 1528  } 

[122]  1529  if (infeasibility) { 

 1530  // Increase estimated degradation to solution 

 1531  estimatedDegradation += CoinMin(object>upEstimate(),object>downEstimate()); 

 1532  numberUnsatisfied_++; 

[460]  1533  sumInfeasibilities_ += infeasibility; 

[122]  1534  // Better priority? Flush choices. 

 1535  if (priorityLevel<bestPriority) { 

 1536  int j; 

 1537  iSmallest=0; 

 1538  for (j=0;j<maximumStrong;j++) { 

 1539  choice[j].upMovement=0.0; 

 1540  delete choice[j].possibleBranch; 

 1541  choice[j].possibleBranch=NULL; 

 1542  } 

 1543  bestPriority = priorityLevel; 

 1544  mostAway=1.0e100; 

 1545  numberStrong=0; 

 1546  } else if (priorityLevel>bestPriority) { 

 1547  continue; 

[44]  1548  } 

[122]  1549  // Check for suitability based on infeasibility. 

 1550  if (infeasibility>mostAway) { 

 1551  //add to list 

 1552  choice[iSmallest].upMovement=infeasibility; 

 1553  delete choice[iSmallest].possibleBranch; 

[640]  1554  CbcSimpleInteger * obj = 

 1555  dynamic_cast <CbcSimpleInteger *>(object) ; 

 1556  if (obj) { 

 1557  choice[iSmallest].possibleBranch=obj>createBranch(solver,&usefulInfo,preferredWay); 

 1558  } else { 

 1559  CbcObject * obj = 

 1560  dynamic_cast <CbcObject *>(object) ; 

 1561  assert (obj); 

 1562  choice[iSmallest].possibleBranch=obj>createBranch(preferredWay); 

 1563  } 

[122]  1564  numberStrong = CoinMax(numberStrong,iSmallest+1); 

 1565  // Save which object it was 

 1566  choice[iSmallest].objectNumber=i; 

 1567  int j; 

 1568  iSmallest=1; 

 1569  mostAway = 1.0e50; 

 1570  for (j=0;j<maximumStrong;j++) { 

 1571  if (choice[j].upMovement<mostAway) { 

 1572  mostAway=choice[j].upMovement; 

 1573  iSmallest=j; 

 1574  } 

[44]  1575  } 

 1576  } 

 1577  } 

 1578  } 

[307]  1579  if (!canDoOneHot&&hotstartSolution) { 

 1580  // switch off as not possible 

 1581  hotstartSolution=NULL; 

 1582  model>setHotstartSolution(NULL,NULL); 

 1583  } 

[122]  1584  if (numberUnsatisfied_) { 

 1585  // some infeasibilities  go to next steps 

[44]  1586  break; 

[122]  1587  } else if (!iPass) { 

 1588  // looks like a solution  get paranoid 

 1589  bool roundAgain=false; 

 1590  // get basis 

 1591  CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver>getWarmStart()); 

 1592  if (!ws) 

[44]  1593  break; 

[122]  1594  for (i=0;i<numberColumns;i++) { 

 1595  double value = saveSolution[i]; 

 1596  if (value<lower[i]) { 

 1597  saveSolution[i]=lower[i]; 

 1598  roundAgain=true; 

 1599  ws>setStructStatus(i,CoinWarmStartBasis::atLowerBound); 

 1600  } else if (value>upper[i]) { 

 1601  saveSolution[i]=upper[i]; 

 1602  roundAgain=true; 

 1603  ws>setStructStatus(i,CoinWarmStartBasis::atUpperBound); 

 1604  } 

[44]  1605  } 

[259]  1606  if (roundAgain&&saveNumberStrong) { 

[122]  1607  // restore basis 

 1608  solver>setWarmStart(ws); 

 1609  delete ws; 

 1610  solver>resolve(); 

 1611  memcpy(saveSolution,solver>getColSolution(),numberColumns*sizeof(double)); 

[135]  1612  model>reserveCurrentSolution(saveSolution); 

[122]  1613  if (!solver>isProvenOptimal()) { 

 1614  // infeasible 

 1615  anyAction=2; 

 1616  break; 

 1617  } 

 1618  } else { 

 1619  delete ws; 

 1620  break; 

 1621  } 

[2]  1622  } 

 1623  } 

[122]  1624  /* Some solvers can do the strong branching calculations faster if 

 1625  they do them all at once. At present only Clp does for ordinary 

 1626  integers but I think this coding would be easy to modify 

 1627  */ 

 1628  bool allNormal=true; // to say if we can do fast strong branching 

 1629  // Say which one will be best 

 1630  int bestChoice=0; 

 1631  double worstInfeasibility=0.0; 

 1632  for (i=0;i<numberStrong;i++) { 

 1633  choice[i].numIntInfeasUp = numberUnsatisfied_; 

 1634  choice[i].numIntInfeasDown = numberUnsatisfied_; 

 1635  choice[i].fix=0; // say not fixed 

 1636  if (!dynamic_cast <const CbcSimpleInteger *> (model>object(choice[i].objectNumber))) 

 1637  allNormal=false; // Something odd so lets skip clever fast branching 

 1638  if ( !model>object(choice[i].objectNumber)>boundBranch()) 

 1639  numberStrong=0; // switch off 

[216]  1640  if ( choice[i].possibleBranch>numberBranches()>2) 

 1641  numberStrong=0; // switch off 

[122]  1642  // Do best choice in case switched off 

 1643  if (choice[i].upMovement>worstInfeasibility) { 

 1644  worstInfeasibility=choice[i].upMovement; 

 1645  bestChoice=i; 

 1646  } 

[6]  1647  } 

[122]  1648  // If we have hit max time don't do strong branching 

 1649  bool hitMaxTime = ( CoinCpuTime()model>getDblParam(CbcModel::CbcStartSeconds) > 

[48]  1650  model>getDblParam(CbcModel::CbcMaximumSeconds)); 

[122]  1651  // also give up if we are looping round too much 

 1652  if (hitMaxTimenumberPassesLeft<=0) 

 1653  numberStrong=0; 

 1654  /* 

 1655  Is strong branching enabled? If so, set up and do it. Otherwise, we'll 

 1656  fall through to simple branching. 

 1657  

 1658  Setup for strong branching involves saving the current basis (for restoration 

 1659  afterwards) and setting up for hot starts. 

 1660  */ 

 1661  if (numberStrong&&saveNumberStrong) { 

 1662  

 1663  bool solveAll=false; // set true to say look at all even if some fixed (experiment) 

 1664  solveAll=true; 

 1665  // worth trying if too many times 

 1666  // Save basis 

 1667  CoinWarmStart * ws = solver>getWarmStart(); 

 1668  // save limit 

 1669  int saveLimit; 

 1670  solver>getIntParam(OsiMaxNumIterationHotStart,saveLimit); 

 1671  if (beforeSolution&&saveLimit<100) 

 1672  solver>setIntParam(OsiMaxNumIterationHotStart,100); // go to end 

[311]  1673  # ifdef COIN_HAS_CLP 

[122]  1674  /* If we are doing all strong branching in one go then we create new arrays 

 1675  to store information. If clp NULL then doing old way. 

 1676  Going down  

 1677  outputSolution[2*i] is final solution. 

 1678  outputStuff[2*i] is status (0  finished, 1 infeas, other unknown 

 1679  outputStuff[2*i+numberStrong] is number iterations 

 1680  On entry newUpper[i] is new upper bound, on exit obj change 

 1681  Going up  

 1682  outputSolution[2*i+1] is final solution. 

 1683  outputStuff[2*i+1] is status (0  finished, 1 infeas, other unknown 

 1684  outputStuff[2*i+1+numberStrong] is number iterations 

[2]  1685  On entry newLower[i] is new lower bound, on exit obj change 

[122]  1686  */ 

 1687  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver); 

 1688  ClpSimplex * clp=NULL; 

 1689  double * newLower = NULL; 

 1690  double * newUpper = NULL; 

 1691  double ** outputSolution=NULL; 

 1692  int * outputStuff=NULL; 

 1693  // Go back to normal way if user wants it 

 1694  if (osiclp&&(osiclp>specialOptions()&16)!=0&&osiclp>specialOptions()>0) 

 1695  allNormal=false; 

[135]  1696  if (osiclp&&!allNormal) { 

 1697  // say do fast 

 1698  int easy=1; 

 1699  osiclp>setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ; 

 1700  } 

[122]  1701  if (osiclp&& allNormal) { 

 1702  clp = osiclp>getModelPtr(); 

 1703  // Clp  do a different way 

 1704  newLower = new double[numberStrong]; 

 1705  newUpper = new double[numberStrong]; 

 1706  outputSolution = new double * [2*numberStrong]; 

 1707  outputStuff = new int [4*numberStrong]; 

 1708  int * which = new int[numberStrong]; 

 1709  int startFinishOptions; 

 1710  int specialOptions = osiclp>specialOptions(); 

 1711  int clpOptions = clp>specialOptions(); 

 1712  int returnCode=0; 

 1713  #define CRUNCH 

 1714  #ifdef CRUNCH 

 1715  // Crunch down problem 

 1716  int numberRows = clp>numberRows(); 

 1717  // Use dual region 

 1718  double * rhs = clp>dualRowSolution(); 

 1719  int * whichRow = new int[3*numberRows]; 

 1720  int * whichColumn = new int[2*numberColumns]; 

 1721  int nBound; 

 1722  ClpSimplex * small = ((ClpSimplexOther *) clp)>crunch(rhs,whichRow,whichColumn,nBound,true); 

 1723  if (!small) { 

 1724  anyAction=2; 

 1725  //printf("XXXX Inf by inspection\n"); 

 1726  delete [] whichColumn; 

 1727  whichColumn=NULL; 

 1728  delete [] whichRow; 

 1729  whichRow=NULL; 

 1730  break; 

 1731  } else { 

 1732  clp = small; 

 1733  } 

 1734  #else 

 1735  int saveLogLevel = clp>logLevel(); 

 1736  int saveMaxIts = clp>maximumIterations(); 

 1737  #endif 

 1738  clp>setLogLevel(0); 

 1739  if((specialOptions&1)==0) { 

 1740  startFinishOptions=0; 

 1741  clp>setSpecialOptions(clpOptions(641024)); 

 1742  } else { 

 1743  startFinishOptions=1+2+4; 

 1744  //startFinishOptions=1+4; // for moment refactorize 

 1745  if((specialOptions&4)==0) 

 1746  clp>setSpecialOptions(clpOptions(6412851210244096)); 

 1747  else 

 1748  clp>setSpecialOptions(clpOptions(64128512102420484096)); 

 1749  } 

 1750  // User may want to clean up before strong branching 

 1751  if ((clp>specialOptions()&32)!=0) { 

 1752  clp>primal(1); 

 1753  if (clp>numberIterations()) 

[640]  1754  model>messageHandler()>message(CBC_ITERATE_STRONG,*model>messagesPointer()) 

[122]  1755  << clp>numberIterations() 

 1756  <<CoinMessageEol; 

 1757  } 

 1758  clp>setMaximumIterations(saveLimit); 

 1759  #ifdef CRUNCH 

 1760  int * backColumn = whichColumn+numberColumns; 

 1761  #endif 

 1762  for (i=0;i<numberStrong;i++) { 

 1763  int iObject = choice[i].objectNumber; 

[640]  1764  const OsiObject * object = model>object(iObject); 

[122]  1765  const CbcSimpleInteger * simple = dynamic_cast <const CbcSimpleInteger *> (object); 

[640]  1766  int iSequence = simple>columnNumber(); 

[122]  1767  newLower[i]= ceil(saveSolution[iSequence]); 

 1768  newUpper[i]= floor(saveSolution[iSequence]); 

 1769  #ifdef CRUNCH 

 1770  iSequence = backColumn[iSequence]; 

 1771  assert (iSequence>=0); 

 1772  #endif 

 1773  which[i]=iSequence; 

 1774  outputSolution[2*i]= new double [numberColumns]; 

 1775  outputSolution[2*i+1]= new double [numberColumns]; 

 1776  } 

 1777  //clp>writeMps("bad"); 

 1778  returnCode=clp>strongBranching(numberStrong,which, 

 1779  newLower, newUpper,outputSolution, 

 1780  outputStuff,outputStuff+2*numberStrong,!solveAll,false, 

 1781  startFinishOptions); 

 1782  #ifndef CRUNCH 

 1783  clp>setSpecialOptions(clpOptions); // restore 

 1784  clp>setMaximumIterations(saveMaxIts); 

 1785  clp>setLogLevel(saveLogLevel); 

 1786  #endif 

 1787  if (returnCode==2) { 

 1788  // bad factorization!!! 

 1789  // Doing normal way 

 1790  // Mark hot start 

 1791  solver>markHotStart(); 

 1792  clp = NULL; 

 1793  } else { 

 1794  #ifdef CRUNCH 

 1795  // extract solution 

 1796  //bool checkSol=true; 

 1797  for (i=0;i<numberStrong;i++) { 

 1798  int iObject = choice[i].objectNumber; 

[640]  1799  const OsiObject * object = model>object(iObject); 

[122]  1800  const CbcSimpleInteger * simple = dynamic_cast <const CbcSimpleInteger *> (object); 

[640]  1801  int iSequence = simple>columnNumber(); 

[122]  1802  which[i]=iSequence; 

 1803  double * sol = outputSolution[2*i]; 

 1804  double * sol2 = outputSolution[2*i+1]; 

 1805  //bool x=true; 

 1806  //bool x2=true; 

 1807  for (int iColumn=numberColumns1;iColumn>=0;iColumn) { 

 1808  int jColumn = backColumn[iColumn]; 

 1809  if (jColumn>=0) { 

 1810  sol[iColumn]=sol[jColumn]; 

 1811  sol2[iColumn]=sol2[jColumn]; 

 1812  } else { 

 1813  sol[iColumn]=saveSolution[iColumn]; 

 1814  sol2[iColumn]=saveSolution[iColumn]; 

 1815  } 

 1816  } 

 1817  } 

 1818  #endif 

 1819  } 

 1820  #ifdef CRUNCH 

 1821  delete [] whichColumn; 

 1822  delete [] whichRow; 

 1823  delete small; 

 1824  #endif 

 1825  delete [] which; 

[32]  1826  } else { 

[122]  1827  // Doing normal way 

 1828  // Mark hot start 

 1829  solver>markHotStart(); 

[32]  1830  } 

[311]  1831  # else /* COIN_HAS_CLP */ 

[277]  1832  

 1833  OsiSolverInterface *clp = NULL ; 

 1834  double **outputSolution = NULL ; 

 1835  int *outputStuff = NULL ; 

 1836  double * newLower = NULL ; 

 1837  double * newUpper = NULL ; 

 1838  

 1839  solver>markHotStart(); 

 1840  

[311]  1841  # endif /* COIN_HAS_CLP */ 

[122]  1842  /* 

 1843  Open a loop to do the strong branching LPs. For each candidate variable, 

 1844  solve an LP with the variable forced down, then up. If a direction turns 

 1845  out to be infeasible or monotonic (i.e., over the dual objective cutoff), 

 1846  force the objective change to be big (1.0e100). If we determine the problem 

 1847  is infeasible, or find a monotone variable, escape the loop. 

 1848  

 1849  TODO: The `restore bounds' part might be better encapsulated as an 

[2]  1850  unbranch() method. Branching objects more exotic than simple integers 

 1851  or cliques might not restrict themselves to variable bounds. 

 1852  

[122]  1853  TODO: Virtuous solvers invalidate the current solution (or give bogus 

[2]  1854  results :) when the bounds are changed out from under them. So we 

 1855  need to do all the work associated with finding a new solution before 

 1856  restoring the bounds. 

[122]  1857  */ 

 1858  for (i = 0 ; i < numberStrong ; i++) 

 1859  { double objectiveChange ; 

 1860  double newObjectiveValue=1.0e100; 

 1861  // status is 0 finished, 1 infeasible and other 

 1862  int iStatus; 

 1863  /* 

 1864  Try the down direction first. (Specify the initial branching alternative as 

 1865  down with a call to way(1). Each subsequent call to branch() performs the 

 1866  specified branch and advances the branch object state to the next branch 

 1867  alternative.) 

 1868  */ 

 1869  if (!clp) { 

 1870  choice[i].possibleBranch>way(1) ; 

 1871  choice[i].possibleBranch>branch() ; 

 1872  bool feasible=true; 

 1873  if (checkFeasibility) { 

 1874  // check branching did not make infeasible 

 1875  int iColumn; 

 1876  int numberColumns = solver>getNumCols(); 

 1877  const double * columnLower = solver>getColLower(); 

 1878  const double * columnUpper = solver>getColUpper(); 

 1879  for (iColumn= 0;iColumn<numberColumns;iColumn++) { 

 1880  if (columnLower[iColumn]>columnUpper[iColumn]+1.0e5) 

 1881  feasible=false; 

 1882  } 

[104]  1883  } 

[122]  1884  if (feasible) { 

 1885  solver>solveFromHotStart() ; 

[259]  1886  numberStrongDone++; 

 1887  numberStrongIterations += solver>getIterationCount(); 

[122]  1888  /* 

 1889  We now have an estimate of objective degradation that we can use for strong 

 1890  branching. If we're over the cutoff, the variable is monotone up. 

 1891  If we actually made it to optimality, check for a solution, and if we have 

 1892  a good one, call setBestSolution to process it. Note that this may reduce the 

 1893  cutoff, so we check again to see if we can declare this variable monotone. 

 1894  */ 

 1895  if (solver>isProvenOptimal()) 

 1896  iStatus=0; // optimal 

 1897  else if (solver>isIterationLimitReached() 

 1898  &&!solver>isDualObjectiveLimitReached()) 

 1899  iStatus=2; // unknown 

 1900  else 

 1901  iStatus=1; // infeasible 

 1902  newObjectiveValue = solver>getObjSense()*solver>getObjValue(); 

 1903  choice[i].numItersDown = solver>getIterationCount(); 

 1904  } else { 

 1905  iStatus=1; // infeasible 

 1906  newObjectiveValue = 1.0e100; 

 1907  choice[i].numItersDown = 0; 

 1908  } 

 1909  } else { 

 1910  iStatus = outputStuff[2*i]; 

 1911  choice[i].numItersDown = outputStuff[2*numberStrong+2*i]; 

[259]  1912  numberStrongDone++; 

 1913  numberStrongIterations += choice[i].numItersDown; 

[122]  1914  newObjectiveValue = objectiveValue+newUpper[i]; 

 1915  solver>setColSolution(outputSolution[2*i]); 

[104]  1916  } 

[268]  1917  objectiveChange = CoinMax(newObjectiveValue  objectiveValue_,0.0); 

[122]  1918  if (!iStatus) { 

 1919  choice[i].finishedDown = true ; 

 1920  if (newObjectiveValue>=model>getCutoff()) { 

 1921  objectiveChange = 1.0e100; // say infeasible 

[259]  1922  numberStrongInfeasible++; 

[122]  1923  } else { 

 1924  // See if integer solution 

 1925  if (model>feasibleSolution(choice[i].numIntInfeasDown, 

[165]  1926  choice[i].numObjInfeasDown) 

 1927  &&model>problemFeasibility()>feasible(model,1)>=0) { 

 1928  model>setBestSolution(CBC_STRONGSOL, 

 1929  newObjectiveValue, 

 1930  solver>getColSolution()) ; 

[264]  1931  // only needed for odd solvers 

 1932  newObjectiveValue = solver>getObjSense()*solver>getObjValue(); 

[268]  1933  objectiveChange = CoinMax(newObjectiveValueobjectiveValue_,0.0) ; 

[181]  1934  model>setLastHeuristic(NULL); 

[137]  1935  model>incrementUsed(solver>getColSolution()); 

[259]  1936  if (newObjectiveValue >= model>getCutoff()) { // *new* cutoff 

[122]  1937  objectiveChange = 1.0e100 ; 

[259]  1938  numberStrongInfeasible++; 

 1939  } 

[165]  1940  } 

[122]  1941  } 

 1942  } else if (iStatus==1) { 

 1943  objectiveChange = 1.0e100 ; 

[259]  1944  numberStrongInfeasible++; 

[122]  1945  } else { 

 1946  // Can't say much as we did not finish 

 1947  choice[i].finishedDown = false ; 

[259]  1948  numberUnfinished++; 

[122]  1949  } 

 1950  choice[i].downMovement = objectiveChange ; 

 1951  

 1952  // restore bounds 

 1953  if (!clp) 

 1954  { for (int j=0;j<numberColumns;j++) { 

 1955  if (saveLower[j] != lower[j]) 

 1956  solver>setColLower(j,saveLower[j]); 

 1957  if (saveUpper[j] != upper[j]) 

 1958  solver>setColUpper(j,saveUpper[j]); 

 1959  } 

 1960  } 

 1961  //printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n", 

 1962  // choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersDown, 

 1963  // choice[i].downMovement,choice[i].finishedDown,choice[i].numIntInfeasDown, 

 1964  // choice[i].numObjInfeasDown); 

 1965  

 1966  // repeat the whole exercise, forcing the variable up 

 1967  if (!clp) { 

 1968  bool feasible=true; 

[216]  1969  // If odd branching then maybe just one possibility 

 1970  if(choice[i].possibleBranch>numberBranchesLeft()>0) { 

 1971  choice[i].possibleBranch>branch(); 

 1972  if (checkFeasibility) { 

 1973  // check branching did not make infeasible 

 1974  int iColumn; 

 1975  int numberColumns = solver>getNumCols(); 

 1976  const double * columnLower = solver>getColLower(); 

 1977  const double * columnUpper = solver>getColUpper(); 

 1978  for (iColumn= 0;iColumn<numberColumns;iColumn++) { 

 1979  if (columnLower[iColumn]>columnUpper[iColumn]+1.0e5) 

 1980  feasible=false; 

 1981  } 

[122]  1982  } 

[216]  1983  } else { 

 1984  // second branch infeasible 

 1985  feasible=false; 

[122]  1986  } 

 1987  if (feasible) { 

 1988  solver>solveFromHotStart() ; 

[259]  1989  numberStrongDone++; 

 1990  numberStrongIterations += solver>getIterationCount(); 

[122]  1991  /* 

 1992  We now have an estimate of objective degradation that we can use for strong 

 1993  branching. If we're over the cutoff, the variable is monotone up. 

 1994  If we actually made it to optimality, check for a solution, and if we have 

 1995  a good one, call setBestSolution to process it. Note that this may reduce the 

 1996  cutoff, so we check again to see if we can declare this variable monotone. 

 1997  */ 

 1998  if (solver>isProvenOptimal()) 

 1999  iStatus=0; // optimal 

 2000  else if (solver>isIterationLimitReached() 

 2001  &&!solver>isDualObjectiveLimitReached()) 

 2002  iStatus=2; // unknown 

 2003  else 

 2004  iStatus=1; // infeasible 

 2005  newObjectiveValue = solver>getObjSense()*solver>getObjValue(); 

 2006  choice[i].numItersUp = solver>getIterationCount(); 

 2007  } else { 

[104]  2008  iStatus=1; // infeasible 

[122]  2009  newObjectiveValue = 1.0e100; 

 2010  choice[i].numItersDown = 0; 

 2011  } 

[104]  2012  } else { 

[122]  2013  iStatus = outputStuff[2*i+1]; 

 2014  choice[i].numItersUp = outputStuff[2*numberStrong+2*i+1]; 

[259]  2015  numberStrongDone++; 

 2016  numberStrongIterations += choice[i].numItersUp; 

[122]  2017  newObjectiveValue = objectiveValue+newLower[i]; 

 2018  solver>setColSolution(outputSolution[2*i+1]); 

[104]  2019  } 

[268]  2020  objectiveChange = CoinMax(newObjectiveValue  objectiveValue_,0.0); 

[122]  2021  if (!iStatus) { 

 2022  choice[i].finishedUp = true ; 

 2023  if (newObjectiveValue>=model>getCutoff()) { 

 2024  objectiveChange = 1.0e100; // say infeasible 

[259]  2025  numberStrongInfeasible++; 

[122]  2026  } else { 

 2027  // See if integer solution 

 2028  if (model>feasibleSolution(choice[i].numIntInfeasUp, 

[165]  2029  choice[i].numObjInfeasUp) 

 2030  &&model>problemFeasibility()>feasible(model,1)>=0) { 

 2031  model>setBestSolution(CBC_STRONGSOL, 

 2032  newObjectiveValue, 

 2033  solver>getColSolution()) ; 

[264]  2034  // only needed for odd solvers 

 2035  newObjectiveValue = solver>getObjSense()*solver>getObjValue(); 

[268]  2036  objectiveChange = CoinMax(newObjectiveValueobjectiveValue_,0.0) ; 

[181]  2037  model>setLastHeuristic(NULL); 

[137]  2038  model>incrementUsed(solver>getColSolution()); 

[259]  2039  if (newObjectiveValue >= model>getCutoff()) { // *new* cutoff 

[122]  2040  objectiveChange = 1.0e100 ; 

[259]  2041  numberStrongInfeasible++; 

 2042  } 

[165]  2043  } 

[122]  2044  } 

 2045  } else if (iStatus==1) { 

 2046  objectiveChange = 1.0e100 ; 

[259]  2047  numberStrongInfeasible++; 

[122]  2048  } else { 

 2049  // Can't say much as we did not finish 

 2050  choice[i].finishedUp = false ; 

[259]  2051  numberUnfinished++; 

[122]  2052  } 

 2053  choice[i].upMovement = objectiveChange ; 

 2054  

 2055  // restore bounds 

 2056  if (!clp) 

 2057  { for (int j=0;j<numberColumns;j++) { 

 2058  if (saveLower[j] != lower[j]) 

 2059  solver>setColLower(j,saveLower[j]); 

 2060  if (saveUpper[j] != upper[j]) 

 2061  solver>setColUpper(j,saveUpper[j]); 

 2062  } 

 2063  } 

 2064  

 2065  //printf("Up on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n", 

 2066  // choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersUp, 

 2067  // choice[i].upMovement,choice[i].finishedUp,choice[i].numIntInfeasUp, 

 2068  // choice[i].numObjInfeasUp); 

 2069  

 2070  /* 

 2071  End of evaluation for this candidate variable. Possibilities are: 

 2072  * Both sides below cutoff; this variable is a candidate for branching. 

 2073  * Both sides infeasible or above the objective cutoff: no further action 

 2074  here. Break from the evaluation loop and assume the node will be purged 

 2075  by the caller. 

 2076  * One side below cutoff: Install the branch (i.e., fix the variable). Break 

 2077  from the evaluation loop and assume the node will be reoptimised by the 

 2078  caller. 

 2079  */ 

[216]  2080  // reset 

 2081  choice[i].possibleBranch>resetNumberBranchesLeft(); 

[122]  2082  if (choice[i].upMovement<1.0e100) { 

 2083  if(choice[i].downMovement<1.0e100) { 

 2084  // feasible  no action 

 2085  } else { 

 2086  // up feasible, down infeasible 

 2087  anyAction=1; 

[135]  2088  //printf("Down infeasible for choice %d sequence %d\n",i, 

 2089  // model>object(choice[i].objectNumber)>columnNumber()); 

[122]  2090  if (!solveAll) { 

 2091  choice[i].possibleBranch>way(1); 

 2092  choice[i].possibleBranch>branch(); 

 2093  break; 

 2094  } else { 

 2095  choice[i].fix=1; 

 2096  } 

 2097  } 

 2098  } else { 

 2099  if(choice[i].downMovement<1.0e100) { 

 2100  // down feasible, up infeasible 

 2101  anyAction=1; 

[135]  2102  //printf("Up infeasible for choice %d sequence %d\n",i, 

 2103  // model>object(choice[i].objectNumber)>columnNumber()); 

[122]  2104  if (!solveAll) { 

 2105  choice[i].possibleBranch>way(1); 

 2106  choice[i].possibleBranch>branch(); 

 2107  break; 

 2108  } else { 

 2109  choice[i].fix=1; 

 2110  } 

 2111  } else { 

 2112  // neither side feasible 

 2113  anyAction=2; 

[135]  2114  //printf("Both infeasible for choice %d sequence %d\n",i, 

 2115  // model>object(choice[i].objectNumber)>columnNumber()); 

[122]  2116  break; 

 2117  } 

 2118  } 

 2119  bool hitMaxTime = ( CoinCpuTime()model>getDblParam(CbcModel::CbcStartSeconds) > 

 2120  model>getDblParam(CbcModel::CbcMaximumSeconds)); 

 2121  if (hitMaxTime) { 

 2122  numberStrong=i+1; 

 2123  break; 

 2124  } 

 2125  } 

 2126  if (!clp) { 

 2127  // Delete the snapshot 

 2128  solver>unmarkHotStart(); 

[2]  2129  } else { 

[122]  2130  delete [] newLower; 

 2131  delete [] newUpper; 

 2132  delete [] outputStuff; 

 2133  int i; 

 2134  for (i=0;i<2*numberStrong;i++) 

 2135  delete [] outputSolution[i]; 

 2136  delete [] outputSolution; 

[2]  2137  } 

[122]  2138  solver>setIntParam(OsiMaxNumIterationHotStart,saveLimit); 

 2139  // restore basis 

 2140  solver>setWarmStart(ws); 

 2141  // Unless infeasible we will carry on 

 2142  // But we could fix anyway 

 2143  if (anyAction==1&&solveAll) { 

 2144  // apply and take off 

 2145  for (i = 0 ; i < numberStrong ; i++) { 

 2146  if (choice[i].fix) { 

 2147  choice[i].possibleBranch>way(choice[i].fix) ; 

 2148  choice[i].possibleBranch>branch() ; 

 2149  } 

 2150  } 

[104]  2151  bool feasible=true; 

 2152  if (checkFeasibility) { 

 2153  // check branching did not make infeasible 

 2154  int iColumn; 

 2155  int numberColumns = solver>getNumCols(); 

 2156  const double * columnLower = solver>getColLower(); 

 2157  const double * columnUpper = solver>getColUpper(); 

 2158  for (iColumn= 0;iColumn<numberColumns;iColumn++) { 

 2159  if (columnLower[iColumn]>columnUpper[iColumn]+1.0e5) 

 2160  feasible=false; 

 2161  } 

 2162  } 

 2163  if (feasible) { 

[122]  2164  // can do quick optimality check 

 2165  int easy=2; 

 2166  solver>setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ; 

 2167  solver>resolve() ; 

 2168  solver>setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ; 

 2169  feasible = solver>isProvenOptimal(); 

 2170  } 

 2171  if (feasible) { 

 2172  memcpy(saveSolution,solver>getColSolution(),numberColumns*sizeof(double)); 

 2173  model>reserveCurrentSolution(saveSolution); 

 2174  memcpy(saveLower,solver>getColLower(),numberColumns*sizeof(double)); 

 2175  memcpy(saveUpper,solver>getColUpper(),numberColumns*sizeof(double)); 

 2176  // Clean up all candidates whih are fixed 

 2177  int numberLeft=0; 

 2178  for (i = 0 ; i < numberStrong ; i++) { 

[135]  2179  CbcStrongInfo thisChoice = choice[i]; 

[122]  2180  choice[i].possibleBranch=NULL; 

[640]  2181  const OsiObject * object = model>object(thisChoice.objectNumber); 

[122]  2182  int preferredWay; 

[640]  2183  double infeasibility = object>infeasibility(&usefulInfo,preferredWay); 

[122]  2184  if (!infeasibility) { 

 2185  // take out 

 2186  delete thisChoice.possibleBranch; 

 2187  } else { 

 2188  choice[numberLeft++]=thisChoice; 

 2189  } 

 2190  } 

 2191  numberStrong=numberLeft; 

 2192  for (;i<maximumStrong;i++) { 

 2193  delete choice[i].possibleBranch; 

 2194  choice[i].possibleBranch=NULL; 

 2195  } 

 2196  // If all fixed then round again 

 2197  if (!numberLeft) { 

 2198  finished=false; 

 2199  numberStrong=0; 

 2200  saveNumberStrong=0; 

 2201  maximumStrong=1; 

 2202  } else { 

 2203  anyAction=0; 

 2204  } 

 2205  // If these two uncommented then different action 

 2206  anyAction=1; 

 2207  finished=true; 

 2208  //printf("some fixed but continuing %d left\n",numberLeft); 

[104]  2209  } else { 

[122]  2210  anyAction=2; // say infeasible 

[104]  2211  } 

[2]  2212  } 

[122]  2213  delete ws; 

[259]  2214  int numberNodes = model>getNodeCount(); 

[262]  2215  // update number of strong iterations etc 

 2216  model>incrementStrongInfo(numberStrongDone,numberStrongIterations, 

 2217  anyAction==2 ? 0:numberStrongInfeasible,anyAction==2); 

[122]  2218  

 2219  /* 

 2220  anyAction >= 0 indicates that strong branching didn't produce any monotone 

 2221  variables. Sift through the candidates for the best one. 

 2222  

 2223  QUERY: Setting numberNodes looks to be a distributed noop. numberNodes is 

 2224  local to this code block. Perhaps should be numberNodes_ from model? 

 2225  Unclear what this calculation is doing. 

 2226  */ 

 2227  if (anyAction>=0) { 

 2228  

 2229  // get average cost per iteration and assume stopped ones 

 2230  // would stop after 50% more iterations at average cost??? !!! ??? 

 2231  double averageCostPerIteration=0.0; 

 2232  double totalNumberIterations=1.0; 

[702]  2233  int smallestNumberInfeasibilities=COIN_INT_MAX; 

[122]  2234  for (i=0;i<numberStrong;i++) { 

 2235  totalNumberIterations += choice[i].numItersDown + 

 2236  choice[i].numItersUp ; 

 2237  averageCostPerIteration += choice[i].downMovement + 

 2238  choice[i].upMovement; 

 2239  smallestNumberInfeasibilities= 

 2240  CoinMin(CoinMin(choice[i].numIntInfeasDown , 

 2241  choice[i].numIntInfeasUp ), 

 2242  smallestNumberInfeasibilities); 

 2243  } 

[197]  2244  //if (smallestNumberInfeasibilities>=numberIntegerInfeasibilities) 

 2245  //numberNodes=1000000; // switch off search for better solution 

[122]  2246  numberNodes=1000000; // switch off anyway 

 2247  averageCostPerIteration /= totalNumberIterations; 

 2248  // all feasible  choose best bet 

 2249  

 2250  // New method does all at once so it can be more sophisticated 

 2251  // in deciding how to balance actions. 

 2252  // But it does need arrays 

 2253  double * changeUp = new double [numberStrong]; 

 2254  int * numberInfeasibilitiesUp = new int [numberStrong]; 

 2255  double * changeDown = new double [numberStrong]; 

 2256  int * numberInfeasibilitiesDown = new int [numberStrong]; 

 2257  CbcBranchingObject ** objects = new CbcBranchingObject * [ numberStrong]; 

 2258  for (i = 0 ; i < numberStrong ; i++) { 

 2259  int iColumn = choice[i].possibleBranch>variable() ; 

[640]  2260  model>messageHandler()>message(CBC_STRONG,*model>messagesPointer()) 

[122]  2261  << i << iColumn 

 2262  <<choice[i].downMovement<<choice[i].numIntInfeasDown 

 2263  <<choice[i].upMovement<<choice[i].numIntInfeasUp 

 2264  <<choice[i].possibleBranch>value() 

 2265  <<CoinMessageEol; 

 2266  changeUp[i]=choice[i].upMovement; 

 2267  numberInfeasibilitiesUp[i] = choice[i].numIntInfeasUp; 

 2268  changeDown[i]=choice[i].downMovement; 

 2269  numberInfeasibilitiesDown[i] = choice[i].numIntInfeasDown; 

 2270  objects[i] = choice[i].possibleBranch; 

 2271  } 

 2272  int whichObject = decision>bestBranch(objects,numberStrong,numberUnsatisfied_, 

 2273  changeUp,numberInfeasibilitiesUp, 

 2274  changeDown,numberInfeasibilitiesDown, 

[222]  2275  objectiveValue_); 

[122]  2276  // move branching object and make sure it will not be deleted 

 2277  if (whichObject>=0) { 

 2278  branch_ = objects[whichObject]; 

[208]  2279  if (model>messageHandler()>logLevel()>3) 

 2280  printf("Choosing column %d\n",choice[whichObject].possibleBranch>variable()) ; 

[122]  2281  choice[whichObject].possibleBranch=NULL; 

 2282  } 

 2283  delete [] changeUp; 

 2284  delete [] numberInfeasibilitiesUp; 

 2285  delete [] changeDown; 

 2286  delete [] numberInfeasibilitiesDown; 

 2287  delete [] objects; 

[2]  2288  } 

[311]  2289  # ifdef COIN_HAS_CLP 

[135]  2290  if (osiclp&&!allNormal) { 

 2291  // back to normal 

 2292  osiclp>setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,NULL) ; 

 2293  } 

[687]  2294  # endif 

[2]  2295  } 

[122]  2296  /* 

 2297  Simple branching. Probably just one, but we may have got here 

 2298  because of an odd branch e.g. a cut 

 2299  */ 

 2300  else { 

 2301  // not strong 

 2302  // C) create branching object 

 2303  branch_ = choice[bestChoice].possibleBranch; 

 2304  choice[bestChoice].possibleBranch=NULL; 

[2]  2305  } 

 2306  } 

 2307  // Set guessed solution value 

 2308  guessedObjectiveValue_ = objectiveValue_+estimatedDegradation; 

 2309  /* 

 2310  Cleanup, then we're outta here. 

 2311  */ 

[640]  2312  if (!model>branchingMethod()dynamicBranchingObject) 

[2]  2313  delete decision; 

 2314  

 2315  for (i=0;i<maximumStrong;i++) 

 2316  delete choice[i].possibleBranch; 

 2317  delete [] choice; 

 2318  delete [] saveLower; 

 2319  delete [] saveUpper; 

 2320  

 2321  // restore solution 

 2322  solver>setColSolution(saveSolution); 

 2323  delete [] saveSolution; 

 2324  return anyAction; 

 2325  } 

 2326  

[135]  2327  /* 

 2328  Version for dynamic pseudo costs. 

[2]  2329  

[135]  2330  **** For now just return if anything odd 

 2331  later allow even if odd 

 2332  

 2333  The routine scans through the object list of the model looking for objects 

 2334  that indicate infeasibility. It tests each object using strong branching 

 2335  and selects the one with the least objective degradation. A corresponding 

 2336  branching object is left attached to lastNode. 

 2337  This version gives preference in evaluation to variables which 

 2338  have not been evaluated many times. It also uses numberStrong 

 2339  to say give up if last few tries have not changed incumbent. 

 2340  See Achterberg, Koch and Martin. 

 2341  

 2342  If strong branching is disabled, a candidate object is chosen essentially 

 2343  at random (whatever object ends up in pos'n 0 of the candidate array). 

 2344  

 2345  If a branching candidate is found to be monotone, bounds are set to fix the 

 2346  variable and the routine immediately returns (the caller is expected to 

 2347  reoptimize). 

 2348  

 2349  If a branching candidate is found to result in infeasibility in both 

 2350  directions, the routine immediately returns an indication of infeasibility. 

 2351  

 2352  Returns: 0 both branch directions are feasible 

 2353  1 branching variable is monotone 

 2354  2 infeasible 

 2355  3 Use another method 

 2356  

 2357  For now just fix on objective from strong branching. 

 2358  */ 

 2359  

[222]  2360  int CbcNode::chooseDynamicBranch (CbcModel *model, CbcNode *lastNode, 

 2361  OsiSolverBranch * & branches,int numberPassesLeft) 

[135]  2362  

 2363  { if (lastNode) 

 2364  depth_ = lastNode>depth_+1; 

 2365  else 

 2366  depth_ = 0; 

[854]  2367  // Go to other choose if hot start 

 2368  if (model>hotstartSolution()) 

 2369  return 3; 

[135]  2370  delete branch_; 

 2371  branch_=NULL; 

 2372  OsiSolverInterface * solver = model>solver(); 

[264]  2373  // get information on solver type 

 2374  const OsiAuxInfo * auxInfo = solver>getAuxiliaryInfo(); 

 2375  const OsiBabSolver * auxiliaryInfo = dynamic_cast<const OsiBabSolver *> (auxInfo); 

[395]  2376  if (!auxiliaryInfo) { 

 2377  // use one from CbcModel 

 2378  auxiliaryInfo = model>solverCharacteristics(); 

 2379  } 

[931]  2380  int numberObjects = model>numberObjects(); 

 2381  bool checkFeasibility = numberObjects>model>numberIntegers(); 

 2382  if ((model>specialOptions()&128)!=0) 

 2383  checkFeasibility=false; // allow 

 2384  // For now return if not simple 

 2385  if (checkFeasibility) 

 2386  return 3; 

[640]  2387  // point to useful information 

 2388  OsiBranchingInformation usefulInfo = model>usefulInformation(); 

 2389  // and modify 

 2390  usefulInfo.depth_=depth_; 

[931]  2391  if ((model>specialOptions()&128)!=0) { 

 2392  // SOS  shadow prices 

 2393  int numberRows = solver>getNumRows(); 

 2394  const double * pi = usefulInfo.pi_; 

 2395  double sumPi=0.0; 

 2396  for (int i=0;i<numberRows;i++) 

 2397  sumPi += fabs(pi[i]); 

 2398  sumPi /= ((double) numberRows); 

 2399  // and scale back 

 2400  sumPi *= 0.01; 

 2401  usefulInfo.defaultDual_ = sumPi; // switch on 

 2402  int numberColumns = solver>getNumCols(); 

 2403  int size = CoinMax(numberColumns,2*numberRows); 

 2404  usefulInfo.usefulRegion_ = new double [size]; 

 2405  CoinZeroN(usefulInfo.usefulRegion_,size); 

 2406  usefulInfo.indexRegion_ = new int [size]; 

 2407  // pi may change 

 2408  usefulInfo.pi_=CoinCopyOfArray(usefulInfo.pi_,numberRows); 

 2409  } 

[395]  2410  assert (auxiliaryInfo); 

[222]  2411  //assert(objectiveValue_ == solver>getObjSense()*solver>getObjValue()); 

[146]  2412  double cutoff =model>getCutoff(); 

[259]  2413  double distanceToCutoff=cutoffobjectiveValue_; 

[135]  2414  const double * lower = solver>getColLower(); 

 2415  const double * upper = solver>getColUpper(); 

[165]  2416  // See what user thinks 

 2417  int anyAction=model>problemFeasibility()>feasible(model,0); 

 2418  if (anyAction) { 

 2419  // will return 2 if infeasible , 0 if treat as integer 

 2420  return anyAction1; 

 2421  } 

[135]  2422  int i; 

[208]  2423  int saveStateOfSearch = model>stateOfSearch(); 

[135]  2424  int numberStrong=model>numberStrong(); 

[264]  2425  if (!auxiliaryInfo>warmStart()) 

 2426  numberStrong=0; 

[135]  2427  // But make more likely to get out after some times 

 2428  int changeStrategy=numberStrong; 

 2429  double changeFactor=1.0; 

 2430  // Use minimum of this and one stored in objects 

 2431  //int numberBeforeTrust = model>numberBeforeTrust(); 

 2432  // Return if doing hot start (in BAB sense) 

[231]  2433  if (model>hotstartSolution()) 

[135]  2434  return 3; 

[833]  2435  //#define RANGING 

[259]  2436  #ifdef RANGING 

[146]  2437  // Pass number 

 2438  int kPass=0; 

[259]  2439  int numberRows = solver>getNumRows(); 

 2440  #endif 

[135]  2441  int numberColumns = model>getNumCols(); 

 2442  double * saveUpper = new double[numberColumns]; 

 2443  double * saveLower = new double[numberColumns]; 

 2444  

 2445  // Save solution in case heuristics need good solution later 

 2446  

 2447  double * saveSolution = new double[numberColumns]; 

 2448  memcpy(saveSolution,solver>getColSolution(),numberColumns*sizeof(double)); 

 2449  model>reserveCurrentSolution(saveSolution); 

 2450  /* 

 2451  Get a branching decision object. Use the default dynamic decision criteria unless 

 2452  the user has loaded a decision method into the model. 

 2453  */ 

 2454  CbcBranchDecision *decision = model>branchingMethod(); 

 2455  if (!decision) 

 2456  decision = new CbcBranchDynamicDecision(); 

[222]  2457  int numberMini=0; 

[210]  2458  int xPen=0; 

 2459  int xMark=0; 

[135]  2460  for (i=0;i<numberColumns;i++) { 

 2461  saveLower[i] = lower[i]; 

 2462  saveUpper[i] = upper[i]; 

 2463  } 

 2464  // Get arrays to sort 

 2465  double * sort = new double[numberObjects]; 

 2466  int * whichObject = new int[numberObjects]; 

[146]  2467  int * objectMark = new int[2*numberObjects+1]; 

[202]  2468  // Arrays with movements 

 2469  double * upEstimate = new double[numberObjects]; 

 2470  double * downEstimate = new double[numberObjects]; 

[135]  2471  CbcStrongInfo * fixObject = new CbcStrongInfo[numberObjects]; 

 2472  double estimatedDegradation=0.0; 

[202]  2473  int numberNodes=model>getNodeCount(); 

[640]  2474  int saveLogLevel = model>logLevel(); 

 2475  if ((numberNodes%500)==0&&false) { 

 2476  model>setLogLevel(6); 

 2477  // Get average up and down costs 

 2478  double averageUp=0.0; 

 2479  double averageDown=0.0; 

 2480  int numberUp=0; 

 2481  int numberDown=0; 

 2482  int i; 

 2483  for ( i=0;i<numberObjects;i++) { 

 2484  OsiObject * object = model>modifiableObject(i); 

 2485  CbcSimpleIntegerDynamicPseudoCost * dynamicObject = 

 2486  dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ; 

 2487  assert(dynamicObject); 

 2488  int numberUp2=0; 

 2489  int numberDown2=0; 

 2490  double up=0.0; 

 2491  double down=0.0; 

 2492  if (dynamicObject>numberTimesUp()) { 

 2493  numberUp++; 

 2494  averageUp += dynamicObject>upDynamicPseudoCost(); 

 2495  numberUp2 += dynamicObject>numberTimesUp(); 

 2496  up = dynamicObject>upDynamicPseudoCost(); 

 2497  } 

 2498  if (dynamicObject>numberTimesDown()) { 

 2499  numberDown++; 

 2500  averageDown += dynamicObject>downDynamicPseudoCost(); 

 2501  numberDown2 += dynamicObject>numberTimesDown(); 

 2502  down = dynamicObject>downDynamicPseudoCost(); 

 2503  } 

 2504  if (numberUp2numberDown2) 

 2505  printf("col %d  up %d times cost %g,  down %d times cost %g\n", 

 2506  dynamicObject>columnNumber(),numberUp2,up,numberDown2,down); 

 2507  } 

 2508  if (numberUp) 

 2509  averageUp /= (double) numberUp; 

 2510  else 

 2511  averageUp=1.0; 

 2512  if (numberDown) 

 2513  averageDown /= (double) numberDown; 

 2514  else 

 2515  averageDown=1.0; 

 2516  printf("total  up %d vars average %g,  down %d vars average %g\n", 

 2517  numberUp,averageUp,numberDown,averageDown); 

 2518  } 

[146]  2519  int numberBeforeTrust = model>numberBeforeTrust(); 

 2520  int numberPenalties = model>numberPenalties(); 

[202]  2521  if (numberBeforeTrust>=1000000) { 

 2522  numberBeforeTrust = numberBeforeTrust % 1000000; 

 2523  numberPenalties=0; 

 2524  } else if (numberBeforeTrust<0) { 

[838]  2525  if (numberBeforeTrust==1) 

 2526  numberPenalties=numberColumns; 

 2527  else if (numberBeforeTrust==2) 

 2528  numberPenalties=0; 

[202]  2529  numberBeforeTrust=0; 

 2530  } 

[135]  2531  // May go round twice if strong branching fixes all local candidates 

 2532  bool finished=false; 

[208]  2533  int numberToFix=0; 

[311]  2534  # ifdef COIN_HAS_CLP 

[208]  2535  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver); 

 2536  int saveClpOptions=0; 

 2537  if (osiclp) { 

 2538  // for faster hot start 

 2539  saveClpOptions = osiclp>specialOptions(); 

[356]  2540  osiclp>setSpecialOptions(saveClpOptions8192); 

[208]  2541  } 

[277]  2542  # else 

 2543  OsiSolverInterface *osiclp = 0 ; 

 2544  # endif 

[640]  2545  const CglTreeProbingInfo * probingInfo = model>probingInfo(); 

[259]  2546  int saveSearchStrategy2 = model>searchStrategy(); 

[931]  2547  #define NODE_NEW 2 

 2548  #ifdef RANGING 

 2549  bool useRanging; 

 2550  #if NODE_NEW !=3 

 2551  useRanging = false; 

 2552  #else 

 2553  useRanging = true; 

 2554  #endif 

 2555  if (saveSearchStrategy2!=2) 

 2556  useRanging=false; 

 2557  #endif 

[259]  2558  if (saveSearchStrategy2<999) { 

[931]  2559  #if NODE_NEW ==4 

 2560  if (saveSearchStrategy2!=2) { 

 2561  #endif 

[259]  2562  // Get average up and down costs 

 2563  double averageUp=0.0; 

 2564  double averageDown=0.0; 

 2565  { 

 2566  int numberUp=0; 

 2567  int numberDown=0; 

 2568  int i; 

 2569  for ( i=0;i<numberObjects;i++) { 

[640]  2570  OsiObject * object = model>modifiableObject(i); 

[259]  2571  CbcSimpleIntegerDynamicPseudoCost * dynamicObject = 

 2572  dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ; 

[931]  2573  if (dynamicObject) { 

 2574  if (dynamicObject>numberTimesUp()) { 

 2575  numberUp++; 

 2576  averageUp += dynamicObject>upDynamicPseudoCost(); 

 2577  } 

 2578  if (dynamicObject>numberTimesDown()) { 

 2579  numberDown++; 

 2580  averageDown += dynamicObject>downDynamicPseudoCost(); 

 2581  } 

 2582  } 

[259]  2583  } 

 2584  if (numberUp) 

 2585  averageUp /= (double) numberUp; 

 2586  else 

 2587  averageUp=1.0; 

 2588  if (numberDown) 

 2589  averageDown /= (double) numberDown; 

 2590  else 

 2591  averageDown=1.0; 

 2592  for ( i=0;i<numberObjects;i++) { 

[640]  2593  OsiObject * object = model>modifiableObject(i); 

[259]  2594  CbcSimpleIntegerDynamicPseudoCost * dynamicObject = 

 2595  dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ; 

[931]  2596  if (dynamicObject) { 

 2597  if (!dynamicObject>numberTimesUp()) 

 2598  dynamicObject>setUpDynamicPseudoCost(averageUp); 

 2599  if (!dynamicObject>numberTimesDown()) 

 2600  dynamicObject>setDownDynamicPseudoCost(averageDown); 

 2601  } 

[259]  2602  } 

 2603  } 

[931]  2604  #if NODE_NEW ==4 

 2605  } else { 

 2606  // pseudo shadow prices 

 2607  model>pseudoShadow(NULL,NULL); 

 2608  } 

 2609  #endif 

[259]  2610  } else if (saveSearchStrategy2<1999) { 

 2611  // pseudo shadow prices 

 2612  model>pseudoShadow(NULL,NULL); 

 2613  } else if (saveSearchStrategy2<2999) { 

 2614  // leave old ones 

 2615  } else if (saveSearchStrategy2<3999) { 

 2616  // pseudo shadow prices at root 

 2617  if (!numberNodes) 

 2618  model>pseudoShadow(NULL,NULL); 

 2619  } else { 

 2620  abort(); 

 2621  } 

 2622  if (saveSearchStrategy2>=0) 

 2623  saveSearchStrategy2 = saveSearchStrategy2 % 1000; 

 2624  if (saveSearchStrategy2==999) 

 2625  saveSearchStrategy2=1; 

 2626  int px[4]={1,1,1,1}; 

 2627  int saveSearchStrategy = saveSearchStrategy2<99 ? saveSearchStrategy2 : saveSearchStrategy2100; 

 2628  bool newWay = saveSearchStrategy2>98; 

 2629  int numberNotTrusted=0; 

[356]  2630  int numberStrongDone=0; 

 2631  int numberUnfinished=0; 

 2632  int numberStrongInfeasible=0; 

 2633  int numberStrongIterations=0; 

[854]  2634  // so we can save lots of news 

 2635  CbcStrongInfo choice; 

 2636  CbcDynamicPseudoCostBranchingObject * choiceObject = NULL; 

 2637  if (model>allDynamic()) { 

 2638  CbcSimpleIntegerDynamicPseudoCost * object = NULL; 

 2639  choiceObject=new CbcDynamicPseudoCostBranchingObject(model,0,1,0.5,object); 

 2640  } 

 2641  choice.possibleBranch=choiceObject; 

[135]  2642  while(!finished) { 

 2643  finished=true; 

 2644  decision>initialize(model); 

 2645  // Some objects may compute an estimate of best solution from here 

 2646  estimatedDegradation=0.0; 

[208]  2647  numberToFix=0; 

[135]  2648  int numberIntegerInfeasibilities=0; // without odd ones 

 2649  int numberToDo=0; 

[140]  2650  int iBestNot=1; 

 2651  int iBestGot=1; 

 2652  double best=0.0; 

[259]  2653  numberNotTrusted=0; 

 2654  numberStrongDone=0; 

 2655  numberUnfinished=0; 

 2656  numberStrongInfeasible=0; 

 2657  numberStrongIterations=0; 

[208]  2658  int * which = objectMark+numberObjects+1; 

 2659  int neededPenalties; 

[640]  2660  int branchingMethod=1; 

[222]  2661  // We may go round this loop three times (only if we think we have solution) 

 2662  for (int iPass=0;iPass<3;iPass++) { 

[135]  2663  

 2664  // compute current state 

 2665  int numberObjectInfeasibilities; // just odd ones 

 2666  model>feasibleSolution( 

 2667  numberIntegerInfeasibilities, 

 2668  numberObjectInfeasibilities); 

 2669  

 2670  // Some objects may compute an estimate of best solution from here 

 2671  estimatedDegradation=0.0; 

 2672  numberUnsatisfied_ = 0; 

[460]  2673  // initialize sum of "infeasibilities" 

 2674  sumInfeasibilities_ = 0.0; 

[702]  2675  int bestPriority=COIN_INT_MAX; 

[640]  2676  int number01 = 0; 

 2677  const fixEntry * entry = NULL; 

 2678  const int * toZero = NULL; 

 2679  const int * toOne = NULL; 

 2680  const int * backward = NULL; 

 2681  int numberUnsatisProbed=0; 

 2682  int numberUnsatisNotProbed=0; // 01 

 2683  if (probingInfo) { 

 2684  number01 = probingInfo>numberIntegers(); 

 2685  entry = probingInfo>fixEntries(); 

 2686  toZero = probingInfo>toZero(); 

 2687  toOne = probingInfo>toOne(); 

 2688  backward = probingInfo>backward(); 

[833]  2689  if (!toZero[number01]number01<numberObjectstrue) { 

[640]  2690  // no info 

 2691  probingInfo=NULL; 

 2692  } 

 2693  } 

[135]  2694  /* 

 2695  Scan for branching objects that indicate infeasibility. Choose candidates 

 2696  using priority as the first criteria, then integer infeasibility. 

 2697  

 2698  The algorithm is to fill the array with a set of good candidates (by 

 2699  infeasibility) with priority bestPriority. Finding a candidate with 

 2700  priority better (less) than bestPriority flushes the choice array. (This 

 2701  serves as initialization when the first candidate is found.) 

 2702  

 2703  */ 

 2704  numberToDo=0; 

[208]  2705  neededPenalties=0; 

[140]  2706  iBestNot=1; 

 2707  double bestNot=0.0; 

 2708  iBestGot=1; 

 2709  best=0.0; 

[819]  2710  /* Problem type as set by user or found by analysis. This will be extended 

 2711  0  not known 

 2712  1  Set partitioning <= 

 2713  2  Set partitioning == 

 2714  3  Set covering 

 2715  4  all + 1 or all +1 and odd 

 2716  */ 

 2717  int problemType = model>problemType(); 

[202]  2718  #define PRINT_STUFF 1 

[135]  2719  for (i=0;i<numberObjects;i++) { 

[640]  2720  OsiObject * object = model>modifiableObject(i); 

[140]  2721  CbcSimpleIntegerDynamicPseudoCost * dynamicObject = 

 2722  dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ; 

[135]  2723  int preferredWay; 

[640]  2724  double infeasibility = object>infeasibility(&usefulInfo,preferredWay); 

[135]  2725  int priorityLevel = object>priority(); 

[259]  2726  #define ZERO_ONE 0 

 2727  #define ZERO_FAKE 1.0e20; 

 2728  #if ZERO_ONE==1 

 2729  // branch on 01 first (temp) 

 2730  if (fabs(saveSolution[dynamicObject>columnNumber()])<1.0) 

 2731  priorityLevel; 

 2732  #endif 

 2733  #if ZERO_ONE==2 

 2734  if (fabs(saveSolution[dynamicObject>columnNumber()])<1.0) 

 2735  infeasibility *= ZERO_FAKE; 

 2736  #endif 

[135]  2737  if (infeasibility) { 

[931]  2738  int iColumn = numberColumns+i; 

 2739  bool gotDown=false; 

 2740  int numberThisDown = 0; 

 2741  bool gotUp=false; 

 2742  int numberThisUp = 0; 

[640]  2743  // check branching method 

[931]  2744  if (dynamicObject) { 

 2745  if (branchingMethod!=dynamicObject>method()) { 

 2746  if (branchingMethod==1) 

 2747  branchingMethod = dynamicObject>method(); 

 2748  else 

 2749  branchingMethod = 100; 

[822]  2750  } 

[931]  2751  iColumn = dynamicObject>columnNumber(); 

 2752  //double gap = saveUpper[iColumn]saveLower[iColumn]; 

 2753  // Give precedence to ones with gap of 1.0 

 2754  //assert(gap>0.0); 

 2755  //infeasibility /= CoinMin(gap,100.0); 

 2756  if (!depth_&&false) { 

 2757  // try closest to 0.5 

 2758  double part =saveSolution[iColumn]floor(saveSolution[iColumn]); 

 2759  infeasibility = fabs(0.5part); 

 2760  } 

 2761  if (problemType>0&&problemType<4&&false) { 

 2762  // try closest to 0.5 

 2763  double part =saveSolution[iColumn]floor(saveSolution[iColumn]); 

 2764  infeasibility = 0.5fabs(0.5part); 

 2765  } 

 2766  if (probingInfo) { 

 2767  int iSeq = backward[iColumn]; 

 2768  assert (iSeq>=0); 

 2769  infeasibility = 1.0 + (toZero[iSeq+1]toZero[iSeq])+ 

 2770  5.0*CoinMin(toOne[iSeq]toZero[iSeq],toZero[iSeq+1]toOne[iSeq]); 

 2771  if (toZero[iSeq+1]>toZero[iSeq]) { 

 2772  numberUnsatisProbed++; 

 2773  } else { 

 2774  numberUnsatisNotProbed++; 

 2775  } 

 2776  } 

 2777  gotDown=false; 

 2778  numberThisDown = dynamicObject>numberTimesDown(); 

 2779  if (numberThisDown>=numberBeforeTrust) 

 2780  gotDown=true; 

 2781  gotUp=false; 

 2782  numberThisUp = dynamicObject>numberTimesUp(); 

 2783  if (numberThisUp>=numberBeforeTrust) 

 2784  gotUp=true; 

 2785  if ((numberNodes%PRINT_STUFF)==0&&PRINT_STUFF>0) 

 2786  printf("%d down %d %g up %d %g  infeas %g\n", 

 2787  i,numberThisDown,object>downEstimate(),numberThisUp,object>upEstimate(), 

 2788  infeasibility); 

 2789  } else { 

 2790  // see if SOS 

 2791  CbcSOS * sosObject = 

 2792  dynamic_cast <CbcSOS *>(object) ; 

 2793  assert (sosObject); 

 2794  gotDown=false; 

 2795  numberThisDown = sosObject>numberTimesDown(); 

 2796  if (numberThisDown>=numberBeforeTrust) 

 2797  gotDown=true; 

 2798  gotUp=false; 

 2799  numberThisUp = sosObject>numberTimesUp(); 

 2800  if (numberThisUp>=numberBeforeTrust) 

 2801  gotUp=true; 

[822]  2802  } 

[135]  2803  // Increase estimated degradation to solution 

 2804  estimatedDegradation += CoinMin(object>upEstimate(),object>downEstimate()); 

[202]  2805  downEstimate[i]=object>downEstimate(); 

 2806  upEstimate[i]=object>upEstimate(); 

[135]  2807  numberUnsatisfied_++; 

[460]  2808  sumInfeasibilities_ += infeasibility; 

[135]  2809  // Better priority? Flush choices. 

 2810  if (priorityLevel<bestPriority) { 

 2811  numberToDo=0; 

 2812  bestPriority = priorityLevel; 

[140]  2813  iBestGot=1; 

 2814  best=0.0; 

[146]  2815  numberNotTrusted=0; 

[135]  2816  } else if (priorityLevel>bestPriority) { 

 2817  continue; 

 2818  } 

[146]  2819  if (!gotUp!gotDown) 

 2820  numberNotTrusted++; 

[135]  2821  // Check for suitability based on infeasibility. 

[146]  2822  if ((gotDown&&gotUp)&&numberStrong>0) { 

[140]  2823  sort[numberToDo]=infeasibility; 

 2824  if (infeasibility>best) { 

 2825  best=infeasibility; 

 2826  iBestGot=numberToDo; 

 2827  } 

 2828  } else { 

[208]  2829  objectMark[neededPenalties]=numberToDo; 

[931]  2830  which[neededPenalties++]=iColumn; 

[259]  2831  sort[numberToDo]=10.0*infeasibility; 

 2832  if (!(numberThisUp+numberThisDown)) 

 2833  sort[numberToDo] *= 100.0; // make even more likely 

[931]  2834  if (iColumn<numberColumns) { 

 2835  double part =saveSolution[iColumn]floor(saveSolution[iColumn]); 

 2836  if (1.0fabs(part0.5)>bestNot) { 

 2837  iBestNot=numberToDo; 

 2838  bestNot = 1.0fabs(part0.5); 

 2839  } 

 2840  } else { 

 2841  // SOS 

 2842  if (sort[numberToDo]>bestNot) { 

 2843  iBestNot=numberToDo; 

 2844  bestNot = sort[numberToDo]; 

 2845  } 

 2846  } 

[140]  2847  } 

[838]  2848  if (model>messageHandler()>logLevel()>3) { 

 2849  printf("%d (%d) down %d %g up %d %g  infeas %g  sort %g solution %g\n", 

 2850  i,iColumn,numberThisDown,object>downEstimate(),numberThisUp,object>upEstimate(), 

 2851  infeasibility,sort[numberToDo],saveSolution[iColumn]); 

 2852  } 

[135]  2853  whichObject[numberToDo++]=i; 

[202]  2854  } else { 

 2855  // for debug 

 2856  downEstimate[i]=1.0; 

 2857  upEstimate[i]=1.0; 

[135]  2858  } 

 2859  } 

 2860  if (numberUnsatisfied_) { 

[640]  2861  if (probingInfo&&false) 

 2862  printf("nunsat %d, %d probed, %d other 01\n",numberUnsatisfied_, 

 2863  numberUnsatisProbed,numberUnsatisNotProbed); 

[135]  2864  // some infeasibilities  go to next steps 

 2865  break; 

 2866  } else if (!iPass) { 

[222]  2867  // may just need resolve 

 2868  solver>resolve(); 

 2869  memcpy(saveSolution,solver>getColSolution(),numberColumns*sizeof(double)); 

 2870  model>reserveCurrentSolution(saveSolution); 

 2871  if (!solver>isProvenOptimal()) { 

 2872  // infeasible 

 2873  anyAction=2; 

 2874  break; 

 2875  } 

 2876  } else if (iPass==1) { 

[135]  2877  // looks like a solution  get paranoid 

 2878  bool roundAgain=false; 

 2879  // get basis 

 2880  CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver>getWarmStart()); 

 2881  if (!ws) 

 2882  break; 

[197]  2883  double tolerance; 

 2884  solver>getDblParam(OsiPrimalTolerance,tolerance); 

[135]  2885  for (i=0;i<numberColumns;i++) { 

 2886  double value = saveSolution[i]; 

[197]  2887  if (value<lower[i]tolerance) { 

[135]  2888  saveSolution[i]=lower[i]; 

 2889  roundAgain=true; 

 2890  ws>setStructStatus(i,CoinWarmStartBasis::atLowerBound); 

[197]  2891  } else if (value>upper[i]+tolerance) { 

[135]  2892  saveSolution[i]=upper[i]; 

 2893  roundAgain=true; 

 2894  ws>setStructStatus(i,CoinWarmStartBasis::atUpperBound); 

 2895  } 

 2896  } 

 2897  if (roundAgain) { 

 2898  // restore basis 

 2899  solver>setWarmStart(ws); 

 2900  solver>setColSolution(saveSolution); 

 2901  delete ws; 

[197]  2902  bool takeHint; 

 2903  OsiHintStrength strength; 

 2904  solver>getHintParam(OsiDoDualInResolve,takeHint,strength); 

 2905  solver>setHintParam(OsiDoDualInResolve,false,OsiHintDo) ; 

[135]  2906  solver>resolve(); 

[197]  2907  solver>setHintParam(OsiDoDualInResolve,takeHint,strength) ; 

[135]  2908  memcpy(saveSolution,solver>getColSolution(),numberColumns*sizeof(double)); 

 2909  model>reserveCurrentSolution(saveSolution); 

 2910  if (!solver>isProvenOptimal()) { 

 2911  // infeasible 

 2912  anyAction=2; 

 2913  break; 

 2914  } 

 2915  } else { 

 2916  delete ws; 

 2917  break; 

 2918  } 

 2919  } 

 2920  } 

 2921  if (anyAction==2) { 

 2922  break; 

 2923  } 

 2924  bool solveAll=false; // set true to say look at all even if some fixed (experiment) 

 2925  solveAll=true; 

[140]  2926  // skip if solution 

 2927  if (!numberUnsatisfied_) 

 2928  break; 

[210]  2929  //bool skipAll = (numberBeforeTrust>20&&numberNodes>20000&&numberNotTrusted==0); 

[687]  2930  bool skipAll = numberNotTrusted==0numberToDo==1; 

[259]  2931  bool doneHotStart=false; 

[640]  2932  int searchStrategy = saveSearchStrategy>=0 ? (saveSearchStrategy%10) : 1; 

[931]  2933  if (0) { 

 2934  int numberIterations = model>getIterationCount(); 

 2935  int numberStrongIterations = model>numberStrongIterations(); 

 2936  printf("node %d saveSearch %d search %d  its %d strong %d\n", 

 2937  numberNodes,saveSearchStrategy,searchStrategy, 

 2938  numberIterations,numberStrongIterations); 

 2939  } 

[259]  2940  #ifndef CBC_WEAK_STRONG 

[640]  2941  if (((numberNodes%20)==0&&searchStrategy!=2)(model>specialOptions()&8)!=0) 

[202]  2942  skipAll=false; 

[259]  2943  #endif 

 2944  if (!newWay) { 

 2945  // 10 up always use %10, 20 up as 10 and allow penalties 

 2946  // But adjust depending on ratio of iterations 

 2947  if (searchStrategy>0&&saveSearchStrategy<10) { 

 2948  if (numberBeforeTrust>=5&&numberBeforeTrust<=10) { 

 2949  if (searchStrategy!=2) { 

[640]  2950  if (depth_>5) { 

 2951  int numberIterations = model>getIterationCount(); 

 2952  int numberStrongIterations = model>numberStrongIterations(); 

 2953  if (numberStrongIterations>numberIterations+10000) { 

 2954  searchStrategy=2; 

 2955  //skipAll=true; 

 2956  } else if (numberStrongIterations*4+1000<numberIterationsdepth_<5) { 

 2957  searchStrategy=3; 

 2958  skipAll=false; 

 2959  } 

 2960  } else { 

 2961  searchStrategy=3; 

 2962  skipAll=false; 

 2963  } 

[259]  2964  } else { 

[640]  2965  //skipAll=true; 

[259]  2966  } 

 2967  } 

 2968  } 

 2969  } else { 

 2970  // But adjust depending on ratio of iterations 

 2971  if (saveSearchStrategy<0) { 

 2972  // unset 

 2973  if ((numberNodes%20)==0(model>specialOptions()&8)!=0) { 

 2974  // Do numberStrong 

 2975  searchStrategy=3; 

 2976  } else if (depth_<5) { 

 2977  // Do numberStrong 

 2978  searchStrategy=2; 

 2979  } else { 

 2980  int numberIterations = model>getIterationCount(); 

 2981  int numberStrongIterations = model>numberStrongIterations(); 

 2982  int numberRows = solver>getNumRows(); 

 2983  if (numberStrongIterations>numberIterations+CoinMin(10000,10*numberRows)) { 

 2984  // off 

 2985  searchStrategy=0; 

 2986  } else if (numberStrongIterations*4+1000<numberIterations) { 

 2987  // Do numberStrong if not trusted 

 2988  searchStrategy=2; 

 2989  } else { 

 2990  searchStrategy=1; 

 2991  } 

 2992  } 

 2993  } 

 2994  if (searchStrategy<3&&(!numberNotTrusted!searchStrategy)) 

 2995  skipAll=true; 

 2996  else 

 2997  skipAll=false; 

 2998  } 

[135]  2999  // worth trying if too many times 

 3000  // Save basis 

[202]  3001  CoinWarmStart * ws = NULL; 

[135]  3002  // save limit 

[202]  3003  int saveLimit=0; 

[259]  3004  solver>getIntParam(OsiMaxNumIterationHotStart,saveLimit); 

[202]  3005  if (!skipAll) { 

 3006  ws = solver>getWarmStart(); 

[259]  3007  int limit=100; 

 3008  #if 0 

 3009  int averageBranchIterations = model>getIterationCount()/(model>getNodeCount()+1); 

 3010  if (numberNodes) 

 3011  limit = CoinMin(CoinMax(limit,2*averageBranchIterations),500); 

 3012  else 

 3013  limit = 500; 

 3014  #endif 

 3015  if ((!saveStateOfSearchsearchStrategy>3)&&saveLimit<limit&&saveLimit==100) 

 3016  solver>setIntParam(OsiMaxNumIterationHotStart,limit); 

[202]  3017  } 

[135]  3018  // Say which one will be best 

 3019  int whichChoice=0; 

[140]  3020  int bestChoice; 

 3021  if (iBestGot>=0) 

 3022  bestChoice=iBestGot; 

 3023  else 

 3024  bestChoice=iBestNot; 

 3025  assert (bestChoice>=0); 

[135]  3026  // If we have hit max time don't do strong branching 

 3027  bool hitMaxTime = ( CoinCpuTime()model>getDblParam(CbcModel::CbcStartSeconds) > 

 3028  model>getDblParam(CbcModel::CbcMaximumSeconds)); 

 3029  // also give up if we are looping round too much 

[640]  3030  if (hitMaxTimenumberPassesLeft<=0(!numberNotTrusted&&false)branchingMethod==11) { 

[140]  3031  int iObject = whichObject[bestChoice]; 

[640]  3032  OsiObject * object = model>modifiableObject(iObject); 

[135]  3033  int preferredWay; 

[640]  3034  object>infeasibility(&usefulInfo,preferredWay); 

 3035  CbcSimpleInteger * obj = 

 3036  dynamic_cast <CbcSimpleInteger *>(object) ; 

 3037  if (obj) { 

 3038  branch_=obj>createBranch(solver,&usefulInfo,preferredWay); 

 3039  } else { 

 3040  CbcObject * obj = 

 3041  dynamic_cast <CbcObject *>(object) ; 

 3042  assert (obj); 

 3043  branch_=obj>createBranch(preferredWay); 

 3044  } 

 3045  { 

 3046  CbcBranchingObject * branchObj = 

 3047  dynamic_cast <CbcBranchingObject *>(branch_) ; 

 3048  assert (branchObj); 

 3049  branchObj>way(preferredWay); 

 3050  } 

[150]  3051  delete ws; 

 3052  ws=NULL; 

[135]  3053  break; 

 3054  } else { 

 3055  // say do fast 

 3056  int easy=1; 

[202]  3057  if (!skipAll) 

 3058  solver>setHintParam(OsiDoInBranchAndCut,true,OsiHintDo,&easy) ; 

[140]  3059  int iDo; 

 3060  #ifdef RANGING 

[931]  3061  if (useRanging) { 

[259]  3062  if ((skipAll&&numberBeforeTrust&&saveSearchStrategy<20)saveSearchStrategy<10) 

[202]  3063  numberPenalties=0; 

[208]  3064  { 

 3065  // off penalties if too much 

 3066  double needed = neededPenalties; 

 3067  needed *= numberRows; 

[259]  3068  if (needed>1.0e6&&numberNodes&&saveSearchStrategy<20) { 

[208]  3069  numberPenalties=0; 

 3070  neededPenalties=0; 

[140]  3071  } 

[208]  3072  } 

[311]  3073  # ifdef COIN_HAS_CLP 

[208]  3074  if (osiclp&&numberPenalties&&neededPenalties) { 

[210]  3075  xPen += neededPenalties; 

[140]  3076  which; 

[208]  3077  which[0]=neededPenalties; 

[140]  3078  osiclp>passInRanges(which); 

 3079  // Mark hot start and get ranges 

[146]  3080  if (kPass) { 

 3081  // until can work out why solution can go funny 

 3082  int save = osiclp>specialOptions(); 

 3083  osiclp>setSpecialOptions(save256); 

 3084  solver>markHotStart(); 

 3085  osiclp>setSpecialOptions(save); 

 3086  } else { 

 3087  solver>markHotStart(); 

 3088  } 

[264]  3089  assert (auxiliaryInfo>warmStart()); 

[259]  3090  doneHotStart=true; 

[210]  3091  xMark++; 

[146]  3092  kPass++; 

[140]  3093  osiclp>passInRanges(NULL); 

[146]  3094  const double * downCost=osiclp>upRange(); 

[150]  3095  const double * upCost=osiclp>downRange(); 

[259]  3096  //printf("numberTodo %d needed %d numberpenalties %d\n",numberToDo,neededPenalties,numberPenalties); 

[208]  3097  double invTrust = 1.0/((double) numberBeforeTrust); 

 3098  for (int i=0;i<neededPenalties;i++) { 

 3099  int j = objectMark[i]; 

 3100  int iObject = whichObject[j]; 

[640]  3101  OsiObject * object = model>modifiableObject(iObject); 

[140]  3102  CbcSimpleIntegerDynamicPseudoCost * dynamicObject = 

 3103  dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ; 

[146]  3104  int iSequence=dynamicObject>columnNumber(); 

[208]  3105  double value = saveSolution[iSequence]; 

 3106  value = floor(value); 

 3107  double upPenalty = CoinMin(upCost[i],1.0e110)*(1.0value); 

 3108  double downPenalty = CoinMin(downCost[i],1.0e110)*value; 

 3109  if (!numberBeforeTrust) { 

 3110  // override 

 3111  downEstimate[iObject]=downPenalty; 

 3112  upEstimate[iObject]=upPenalty; 

 3113  } else { 

 3114  int numberThisDown = dynamicObject>numberTimesDown(); 

 3115  if (numberThisDown<numberBeforeTrust) { 

 3116  double fraction = ((double) numberThisDown)*invTrust; 

 3117  downEstimate[iObject] = fraction*downEstimate[iObject]+(1.0fraction)*downPenalty; 

[202]  3118  } 

[208]  3119  int numberThisUp = dynamicObject>numberTimesUp(); 

 3120  if (numberThisUp<numberBeforeTrust) { 

 3121  double fraction = ((double) numberThisUp)*invTrust; 

 3122  upEstimate[iObject] = fraction*upEstimate[iObject]+(1.0fraction)*upPenalty; 

[140]  3123  } 

 3124  } 

[208]  3125  sort[j] =  CoinMin(downEstimate[iObject],upEstimate[iObject]); 

[259]  3126  #ifdef CBC_WEAK_STRONG 

 3127  sort[j] = 1.0e10; // make more likely to be chosen 

 3128  #endif 

 3129  //if ((numberNodes%PRINT_STUFF)==0&&PRINT_STUFF>0) 

 3130  if (!numberNodes) 

[208]  3131  printf("%d pen down ps %g > %g up ps %g > %g\n", 

 3132  iObject,downCost[i],downPenalty,upCost[i],upPenalty); 

[140]  3133  } 

[277]  3134  } else 

[311]  3135  # endif /* COIN_HAS_CLP */ 

[277]  3136  { 

[202]  3137  if (!skipAll) { 

 3138  // Mark hot start 

 3139  solver>markHotStart(); 

[259]  3140  doneHotStart=true; 

[264]  3141  assert (auxiliaryInfo>warmStart()); 

[210]  3142  xMark++; 

[208]  3143  //if (solver>isProvenPrimalInfeasible()) 

 3144  //printf("**** Hot start says node infeasible\n"); 

[202]  3145  } 

[140]  3146  // make sure best will be first 

 3147  if (iBestGot>=0) 

[202]  3148  sort[iBestGot]=1.0e120; 

[140]  3149  } 

[931]  3150  } else { 

 3151  #endif /* RANGING */ 

[259]  3152  if (!skipAll) { 

 3153  // Mark hot start 

 3154  doneHotStart=true; 

[264]  3155  assert (auxiliaryInfo>warmStart()); 

[259]  3156  solver>markHotStart(); 

 3157  xMark++; 

 3158  } 

[140]  3159  // make sure best will be first 

 3160  if (iBestGot>=0) 

 3161  sort[iBestGot]=COIN_DBL_MAX; 

[931]  3162  #ifdef RANGING 

 3163  } 

[277]  3164  #endif /* RANGING */ 

[146]  3165  // Actions 0  exit for repeat, 1 resolve and try old choice,2 exit for continue 

 3166  #define ACTION 0 

 3167  #if ACTION<2 

 3168  if (anyAction) 

 3169  numberToDo=0; // skip as we will be trying again 

 3170  #endif 

[208]  3171  // Sort 

 3172  CoinSort_2(sort,sort+numberToDo,whichObject); 

 3173  // Change in objective opposite infeasible 

 3174  double worstFeasible=0.0; 

[146]  3175  // Just first if strong off 

 3176  if (!numberStrong) 

 3177  numberToDo=CoinMin(numberToDo,1); 

[931]  3178  #if NODE_NEW 

 3179  if (searchStrategy==2) 

 3180  numberToDo=CoinMin(numberToDo,20); 

 3181  #endif 

[140]  3182  iDo=0; 

 3183  int saveLimit2; 

 3184  solver>getIntParam(OsiMaxNumIterationHotStart,saveLimit2); 

[202]  3185  bool doQuickly = false; // numberToDo>2*numberStrong; 

[640]  3186  if (searchStrategy==2) 

 3187  doQuickly=true; 

[201]  3188  //printf("todo %d, strong %d\n",numberToDo,numberStrong); 

 3189  int numberTest=numberNotTrusted>0 ? numberStrong : (numberStrong+1)/2; 

[259]  3190  int numberTest2 = 2*numberStrong; 

[838]  3191  //double distanceToCutoff2 = model>getCutoff()objectiveValue_; 

[259]  3192  if (!newWay) { 

 3193  if (searchStrategy==3) { 

 3194  // Previously decided we need strong 

 3195  doQuickly=false; 

 3196  numberTest = numberStrong; 

 3197  //numberTest2 = 1000000; 

 3198  } 

[708]  3199  //if (searchStrategy<0searchStrategy==1) 

[259]  3200  //numberTest2 = 1000000; 

[202]  3201  #if 0 

 3202  if (numberBeforeTrust>20&&(numberNodes>20000(numberNodes>200&&numberNotTrusted==0))) { 

 3203  if ((numberNodes%20)!=0) { 

 3204  numberTest=0; 

 3205  doQuickly=true; 

 3206  } 

 3207  } 

 3208  #else 

 3209  // Try nearly always off 

[259]  3210  if (searchStrategy<2) { 

 3211  if ((numberNodes%20)!=0) { 

 3212  if ((model>specialOptions()&8)==0) { 

 3213  numberTest=0; 

 3214  doQuickly=true; 

 3215  } 

 3216  } else { 

 3217  doQuickly=false; 

 3218  numberTest=2*numberStrong; 

[356]  3219  skipAll=false; 

[259]  3220  } 

 3221  } else if (searchStrategy!=3) { 

[202]  3222  doQuickly=true; 

[259]  3223  numberTest=numberStrong; 

[202]  3224  } 

 3225  #endif 

[640]  3226  if (depth_<8&&numberStrong) { 

[259]  3227  if (searchStrategy!=2) { 

 3228  doQuickly=false; 

[640]  3229  int numberRows = solver>getNumRows(); 

 3230  // whether to do this or not is important  think 

 3231  if (numberRows<300numberRows+numberColumns<2500) { 

 3232  if (depth_<7) 

 3233  numberStrong = CoinMin(3*numberStrong,numberToDo); 

 3234  if (!depth_) 

 3235  numberStrong=CoinMin(6*numberStrong,numberToDo); 

 3236  } 

[259]  3237  numberTest=numberStrong; 

[356]  3238  skipAll=false; 

[259]  3239  } 

[833]  3240  //model>setStateOfSearch(2); // use min min 

[208]  3241  } 

 3242  // could adjust using average iterations per branch 

 3243  // double average = ((double)model>getIterationCount())/ 

 3244  //((double) model>getNodeCount()+1.0); 

[202]  3245  // if too many and big then just do 10 its 

[208]  3246  if (!skipAll&&saveStateOfSearch) { 

[259]  3247  //if (numberNotTrusted>3*numberStrong&&numberRows>250&&numberColumns>1000&&saveLimit==100) 

[208]  3248  // off solver>setIntParam(OsiMaxNumIterationHotStart,10); 

[202]  3249  } 

[208]  3250  // make negative for test 

 3251  distanceToCutoff =  distanceToCutoff; 

 3252  if (numberObjects>100) { 

 3253  // larger 

 3254  distanceToCutoff *= 100.0; 

 3255  } 

[202]  3256  distanceToCutoff = COIN_DBL_MAX; 

[208]  3257  // Do at least 5 strong 

 3258  if (numberColumns<1000&&(depth_<15numberNodes<1000000)) 

 3259  numberTest = CoinMax(numberTest,5); 

[210]  3260  if ((model>specialOptions()&8)==0) { 

 3261  if (skipAll) { 

 3262  numberTest=0; 

 3263  doQuickly=true; 

 3264  } 

 3265  } else { 

 3266  // do 5 as strong is fixing 

 3267  numberTest = CoinMax(numberTest,5); 

 3268  } 

[259]  3269  } else { 

[931]  3270  if (!depth_&&numberToDo<200&&solver>getNumRows()<2000&& 

 3271  numberColumns<2000&&false) { 

 3272  numberStrong = numberToDo; 

 3273  doQuickly = false; 

 3274  } 

[259]  3275  int numberTest=numberNotTrusted>0 ? numberStrong : (numberStrong+1)/2; 

 3276  int numberTest2 = 2*numberStrong; 

 3277  if (searchStrategy>=3) { 

 3278  // Previously decided we need strong 

 3279  doQuickly=false; 

 3280  if (depth_<7) 

 3281  numberStrong *=3; 

 3282  if (!depth_) 

 3283  numberStrong=CoinMin(6*numberStrong,numberToDo); 

 3284  numberTest = numberStrong; 

 3285  numberTest2 *= 2; 

 3286  } else if (searchStrategy==2(searchStrategy==1&&depth_<6)) { 

 3287  numberStrong *=2; 

 3288  if (!depth_) 

 3289  numberStrong=CoinMin(2*numberStrong,numberToDo); 

 3290  numberTest = numberStrong; 

 3291  } else if (searchStrategy==1&&numberNotTrusted) { 

 3292  numberTest = numberStrong; 

 3293  } else { 

 3294  numberTest=0; 

 3295  skipAll=true; 

 3296  } 

 3297  distanceToCutoff=model>getCutoff()objectiveValue_; 

 3298  // make negative for test 

 3299  distanceToCutoff =  distanceToCutoff; 

 3300  if (numberObjects>100) { 

 3301  // larger 

 3302  distanceToCutoff *= 100.0; 

 3303  } 

 3304  distanceToCutoff = COIN_DBL_MAX; 

 3305  if (skipAll) { 

 3306  numberTest=0; 

 3307  doQuickly=true; 

 3308  } 

 3309  } 

[356]  3310  // temp  always switch off 

 3311  if (0) { 

 3312  int numberIterations = model>getIterationCount(); 

 3313  int numberStrongIterations = model>numberStrongIterations(); 

[931]  3314  int numberRows = solver>getNumRows(); 

 3315  if (numberStrongIterations>numberIterations+CoinMin(100,10*numberRows)&&depth_>=5) { 

[356]  3316  skipAll=true; 

 3317  newWay=false; 

 3318  numberTest=0; 

 3319  doQuickly=true; 

 3320  } 

 3321  } 

[931]  3322  #if 0 

[356]  3323  // temp  always switch on 

 3324  if (0) { 

 3325  int numberIterations = model>getIterationCount(); 

 3326  int numberStrongIterations = model>numberStrongIterations(); 

 3327  if (2*numberStrongIterations<numberIterationsdepth_<=5) { 

 3328  skipAll=false; 

 3329  newWay=false; 

 3330  numberTest=CoinMax(numberTest,numberStrong); 

 3331  doQuickly=false; 

 3332  } 

 3333  } 

 3334  #endif 

[259]  3335  px[0]=numberTest; 

 3336  px[1]=numberTest2; 

 3337  px[2]= doQuickly ? 1 : 1; 

 3338  px[3]=numberStrong; 

[640]  3339  if (!newWay) { 

 3340  if (numberColumns>8*solver>getNumRows()&&false) { 

 3341  printf("skipAll %c doQuickly %c numberTest %d numberTest2 %d numberNot %d\n", 

 3342  skipAll ? 'Y' : 'N',doQuickly ? 'Y' : 'N',numberTest,numberTest2,numberNotTrusted); 

 3343  numberTest = CoinMin(numberTest,model>numberStrong()); 

 3344  numberTest2 = CoinMin(numberTest2,model>numberStrong()); 

 3345  printf("new test,test2 %d %d\n",numberTest,numberTest2); 

 3346  } 

 3347  } 

[259]  3348  //printf("skipAll %c doQuickly %c numberTest %d numberTest2 %d numberNot %d\n", 

[640]  3349  // skipAll ? 'Y' : 'N',doQuickly ? 'Y' : 'N',numberTest,numberTest2,numberNotTrusted); 

[222]  3350  // See if we want mini tree 

 3351  bool wantMiniTree=false; 

 3352  if (model>sizeMiniTree()&&depth_>7&&saveStateOfSearch>0) 

 3353  wantMiniTree=true; 

 3354  numberMini=0; 

[356]  3355  //if (skipAll&&numberTest==0&&doQuickly) 

 3356  //numberToDo = 1; // trust previous stuff 

 3357  bool couldChooseFirst = false ; //(skipAll&&numberTest==0&&doQuickly); 

[640]  3358  //skipAll=false; 

[931]  3359  int realMaxHotIterations=999999; 

 3360  #if 0 

 3361  if (saveSearchStrategy==1&&!model>parentModel()&& 

 3362  !depth_&&saveLimit==100) { 

 3363  realMaxHotIterations=saveLimit; 

 3364  saveLimit2=200; 

 3365  solver>setIntParam(OsiMaxNumIterationHotStart,saveLimit2); 

 3366  } 

 3367  #endif 

[135]  3368  for ( iDo=0;iDo<numberToDo;iDo++) { 

 3369  int iObject = whichObject[iDo]; 

[640]  3370  OsiObject * object = model>modifiableObject(iObject); 

[135]  3371  CbcSimpleIntegerDynamicPseudoCost * dynamicObject = 

 3372  dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ; 

[931]  3373  int iColumn = dynamicObject ? dynamicObject>columnNumber() : numberColumns+iObject; 

[135]  3374  int preferredWay; 

[640]  3375  double infeasibility = object>infeasibility(&usefulInfo,preferredWay); 

 3376  // may have become feasible 

 3377  if (!infeasibility) 

 3378  continue; 

 3379  CbcSimpleInteger * obj = 

 3380  dynamic_cast <CbcSimpleInteger *>(object) ; 

 3381  if (obj) { 

[854]  3382  if (choiceObject) { 

 3383  obj>fillCreateBranch(choiceObject,&usefulInfo,preferredWay); 

 3384  choiceObject>setObject(dynamicObject); 

 3385  } else { 

 3386  choice.possibleBranch=obj>createBranch(solver,&usefulInfo,preferredWay); 

 3387  } 

[640]  3388  } else { 

 3389  CbcObject * obj = 

 3390  dynamic_cast <CbcObject *>(object) ; 

 3391  assert (obj); 

 3392  choice.possibleBranch=obj>createBranch(preferredWay); 

 3393  } 

[135]  3394  // Save which object it was 

 3395  choice.objectNumber=iObject; 

 3396  choice.numIntInfeasUp = numberUnsatisfied_; 

 3397  choice.numIntInfeasDown = numberUnsatisfied_; 

[202]  3398  choice.upMovement = upEstimate[iObject]; 

 3399  choice.downMovement = downEstimate[iObject]; 

 3400  assert (choice.upMovement>=0.0); 

 3401  assert (choice.downMovement>=0.0); 

[135]  3402  choice.fix=0; // say not fixed 

[640]  3403  double maxChange = 0.5*(choice.upMovement+choice.downMovement); 

 3404  maxChange = CoinMin(choice.upMovement,choice.downMovement); 

 3405  maxChange = CoinMax(choice.upMovement,choice.downMovement); 

 3406  if (searchStrategy==2) 

 3407  maxChange = COIN_DBL_MAX; 

 3408  //maxChange *= 5.0; 

[931]  3409  if (dynamicObject&&dynamicObject>method()==1) 

[640]  3410  maxChange *= 0.1; // probing 

[135]  3411  // see if can skip strong branching 

 3412  int canSkip = choice.possibleBranch>fillStrongInfo(choice); 

[833]  3413  #if 0 

[259]  3414  if (!newWay) { 

[640]  3415  if ((maxChange>distanceToCutoff2)&&(!doQuickly(numberTest>0&&searchStrategy!=2))) 

[200]  3416  canSkip=0; 

[259]  3417  } else { 

 3418  if (skipAll) 

 3419  canSkip=1; 

 3420  else if (numberTest>0&&searchStrategy>=3) 

 3421  canSkip=0; 

 3422  } 

[202]  3423  if (!numberBeforeTrust) { 

 3424  canSkip=1; 

 3425  } 

 3426  if (sort[iDo]<distanceToCutoff) 

 3427  canSkip=0; 

[259]  3428  if (((numberTest2<=0&&numberTest<=0)skipAll)&&sort[iDo]>distanceToCutoff) { 

[202]  3429  canSkip=1; // always skip 

[293]  3430  if (iDo>20) { 

[854]  3431  if (!choiceObject) { 

 3432  delete choice.possibleBranch; 

 3433  choice.possibleBranch=NULL; 

 3434  } 

[259]  3435  break; // give up anyway 

[293]  3436  } 

[259]  3437  } 

[833]  3438  #else 

 3439  if (((numberTest2<=0&&numberTest<=0)skipAll)&&sort[iDo]>distanceToCutoff) { 

 3440  //canSkip=1; // always skip 

 3441  if (iDo>20) { 

[854]  3442  if (!choiceObject) { 

 3443  delete choice.possibleBranch; 

 3444  choice.possibleBranch=NULL; 

 3445  } 

[833]  3446  break; // give up anyway 

 3447  } 

 3448  } 

 3449  #endif 

[931]  3450  if (model>messageHandler()>logLevel()>3&&numberBeforeTrust&&dynamicObject) 

[135]  3451  dynamicObject>print(1,choice.possibleBranch>value()); 

[140]  3452  // was if (!canSkip) 

[259]  3453  if (newWay) 

 3454  numberTest2; 

[135]  3455  if (!canSkip) { 

[356]  3456  //#ifndef RANGING 

 3457  if (!doneHotStart) { 

 3458  // Mark hot start 

 3459  doneHotStart=true; 

 3460  assert (auxiliaryInfo>warmStart()); 

 3461  solver>markHotStart(); 

 3462  xMark++; 

 3463  } 

 3464  //#endif 

 3465  assert (!couldChooseFirst); 

[201]  3466  numberTest; 

[259]  3467  if (!newWay) 

[202]  3468  numberTest2; 

[140]  3469  // just do a few 

[931]  3470  #if NODE_NEW == 5  NODE_NEW == 2 

[259]  3471  //if (canSkip) 

[931]  3472  if (searchStrategy==2) 

 3473  solver>setIntParam(OsiMaxNumIterationHotStart,10); 

 3474  #endif 

[135]  3475  double objectiveChange ; 

 3476  double newObjectiveValue=1.0e100; 

 3477  int j; 

 3478  // status is 0 finished, 1 infeasible and other 

 3479  int iStatus; 

[687]  3480  if (0) { 

 3481  CbcDynamicPseudoCostBranchingObject * cbcobj = dynamic_cast<CbcDynamicPseudoCostBranchingObject *> (choice.possibleBranch); 

 3482  if (cbcobj) { 

 3483  CbcSimpleIntegerDynamicPseudoCost * object = cbcobj>object(); 

 3484  printf("strong %d ",iDo); 

 3485  object>print(1,0.5); 

 3486  } 

 3487  } 

[135]  3488  /* 

 3489  Try the down direction first. (Specify the initial branching alternative as 

 3490  down with a call to way(1). Each subsequent call to branch() performs the 

 3491  specified branch and advances the branch object state to the next branch 

 3492  alternative.) 

 3493  */ 

 3494  choice.possibleBranch>way(1) ; 

[687]  3495  #if NEW_UPDATE_OBJECT==0 

[135]  3496  decision>saveBranchingObject( choice.possibleBranch); 

[687]  3497  #endif 

[135]  3498  choice.possibleBranch>branch() ; 

 3499  solver>solveFromHotStart() ; 

[264]  3500  bool needHotStartUpdate=false; 

[259]  3501  numberStrongDone++; 

 3502  numberStrongIterations += solver>getIterationCount(); 

[135]  3503  /* 

 3504  We now have an estimate of objective degradation that we can use for strong 

 3505  branching. If we're over the cutoff, the variable is monotone up. 

 3506  If we actually made it to optimality, check for a solution, and if we have 

 3507  a good one, call setBestSolution to process it. Note that this may reduce the 

 3508  cutoff, so we check again to see if we can declare this variable monotone. 

 3509  */ 

 3510  if (solver>isProvenOptimal()) 

 3511  iStatus=0; // optimal 

 3512  else if (solver>isIterationLimitReached() 

 3513  &&!solver>isDualObjectiveLimitReached()) 

 3514  iStatus=2; // unknown 

 3515  else 

 3516  iStatus=1; // infeasible 

[931]  3517  if (iStatus!=2&&solver>getIterationCount()> 

 3518  realMaxHotIterations) 

 3519  numberUnfinished++; 

[135]  3520  newObjectiveValue = solver>getObjSense()*solver>getObjValue(); 

 3521  choice.numItersDown = solver>getIterationCount(); 

[268]  3522  objectiveChange = CoinMax(newObjectiveValue  objectiveValue_,0.0); 

[687]  3523  // Update branching information if wanted 

 3524  #if NEW_UPDATE_OBJECT==0 

[135]  3525  decision>updateInformation( solver,this); 

[687]  3526  #elif NEW_UPDATE_OBJECT<2 

 3527  CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch); 

 3528  if (cbcobj) { 

 3529  CbcObject * object = cbcobj>object(); 

[931]  3530  assert (object); 

[687]  3531  CbcObjectUpdateData update = object>createUpdateInformation(solver,this,cbcobj); 

 3532  object>updateInformation(update); 

 3533  } else { 

 3534  decision>updateInformation( solver,this); 

 3535  } 

 3536  #else 

 3537  CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch); 

 3538  if (cbcobj) { 

 3539  CbcObject * object = cbcobj>object(); 

[931]  3540  assert (object) ; 

[687]  3541  CbcObjectUpdateData update = object>createUpdateInformation(solver,this,cbcobj); 

 3542  update.objectNumber_ = choice.objectNumber; 

 3543  model>addUpdateInformation(update); 

 3544  } else { 

 3545  decision>updateInformation( solver,this); 

