Changeset 2469 for trunk/Cbc/examples/CbcBranchLink.cpp
 Timestamp:
 Jan 6, 2019 6:17:46 PM (10 months ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/Cbc/examples/CbcBranchLink.cpp
r1900 r2469 17 17 #include "CoinPackedMatrix.hpp" 18 18 19 // Default Constructor 20 CbcLink::CbcLink 21 : CbcObject() ,22 weights_(NULL),23 numberMembers_(0),24 numberLinks_(0),25 which_(NULL),26 19 // Default Constructor 20 CbcLink::CbcLink() 21 : CbcObject() 22 , weights_(NULL) 23 , numberMembers_(0) 24 , numberLinks_(0) 25 , which_(NULL) 26 , sosType_(1) 27 27 { 28 28 } 29 29 30 30 // Useful constructor (which are indices) 31 CbcLink::CbcLink (CbcModel * model,int numberMembers,32 int numberLinks, int first , const double *weights, int identifier)33 : CbcObject(model) ,34 numberMembers_(numberMembers),35 numberLinks_(numberLinks),36 which_(NULL),37 38 { 39 id_ =identifier;31 CbcLink::CbcLink(CbcModel *model, int numberMembers, 32 int numberLinks, int first, const double *weights, int identifier) 33 : CbcObject(model) 34 , numberMembers_(numberMembers) 35 , numberLinks_(numberLinks) 36 , which_(NULL) 37 , sosType_(1) 38 { 39 id_ = identifier; 40 40 if (numberMembers_) { 41 41 weights_ = new double[numberMembers_]; 42 which_ = new int[numberMembers_ *numberLinks_];42 which_ = new int[numberMembers_ * numberLinks_]; 43 43 if (weights) { 44 memcpy(weights_, weights,numberMembers_*sizeof(double));44 memcpy(weights_, weights, numberMembers_ * sizeof(double)); 45 45 } else { 46 for (int i =0;i<numberMembers_;i++)47 weights_[i] =i;46 for (int i = 0; i < numberMembers_; i++) 47 weights_[i] = i; 48 48 } 49 49 // weights must be increasing 50 50 int i; 51 for (i =0;i<numberMembers_;i++)52 assert (i == 0  weights_[i]>weights_[i1]+1.0e12);53 for (i =0;i<numberMembers_*numberLinks_;i++) {54 which_[i] =first+i;51 for (i = 0; i < numberMembers_; i++) 52 assert(i == 0  weights_[i] > weights_[i  1] + 1.0e12); 53 for (i = 0; i < numberMembers_ * numberLinks_; i++) { 54 which_[i] = first + i; 55 55 } 56 56 } else { … … 60 60 61 61 // Useful constructor (which are indices) 62 CbcLink::CbcLink (CbcModel * model,int numberMembers,63 int numberLinks, int sosType, const int * which , const double *weights, int identifier)64 : CbcObject(model) ,65 numberMembers_(numberMembers),66 numberLinks_(numberLinks),67 which_(NULL),68 69 { 70 id_ =identifier;62 CbcLink::CbcLink(CbcModel *model, int numberMembers, 63 int numberLinks, int sosType, const int *which, const double *weights, int identifier) 64 : CbcObject(model) 65 , numberMembers_(numberMembers) 66 , numberLinks_(numberLinks) 67 , which_(NULL) 68 , sosType_(sosType) 69 { 70 id_ = identifier; 71 71 if (numberMembers_) { 72 72 weights_ = new double[numberMembers_]; 73 which_ = new int[numberMembers_ *numberLinks_];73 which_ = new int[numberMembers_ * numberLinks_]; 74 74 if (weights) { 75 memcpy(weights_, weights,numberMembers_*sizeof(double));75 memcpy(weights_, weights, numberMembers_ * sizeof(double)); 76 76 } else { 77 for (int i =0;i<numberMembers_;i++)78 weights_[i] =i;77 for (int i = 0; i < numberMembers_; i++) 78 weights_[i] = i; 79 79 } 80 80 // weights must be increasing 81 81 int i; 82 for (i =0;i<numberMembers_;i++)83 assert (i == 0  weights_[i]>weights_[i1]+1.0e12);84 for (i =0;i<numberMembers_*numberLinks_;i++) {85 which_[i] = which[i];82 for (i = 0; i < numberMembers_; i++) 83 assert(i == 0  weights_[i] > weights_[i  1] + 1.0e12); 84 for (i = 0; i < numberMembers_ * numberLinks_; i++) { 85 which_[i] = which[i]; 86 86 } 87 87 } else { … … 90 90 } 91 91 92 // Copy constructor 93 CbcLink::CbcLink ( const CbcLink &rhs)94 : CbcObject(rhs)92 // Copy constructor 93 CbcLink::CbcLink(const CbcLink &rhs) 94 : CbcObject(rhs) 95 95 { 96 96 numberMembers_ = rhs.numberMembers_; … … 98 98 sosType_ = rhs.sosType_; 99 99 if (numberMembers_) { 100 weights_ = CoinCopyOfArray(rhs.weights_, numberMembers_);101 which_ = CoinCopyOfArray(rhs.which_, numberMembers_*numberLinks_);100 weights_ = CoinCopyOfArray(rhs.weights_, numberMembers_); 101 which_ = CoinCopyOfArray(rhs.which_, numberMembers_ * numberLinks_); 102 102 } else { 103 103 weights_ = NULL; … … 113 113 } 114 114 115 // Assignment operator 116 CbcLink & 117 CbcLink::operator=( const CbcLink&rhs)118 { 119 if (this !=&rhs) {115 // Assignment operator 116 CbcLink & 117 CbcLink::operator=(const CbcLink &rhs) 118 { 119 if (this != &rhs) { 120 120 CbcObject::operator=(rhs); 121 delete 122 delete 121 delete[] weights_; 122 delete[] which_; 123 123 numberMembers_ = rhs.numberMembers_; 124 124 numberLinks_ = rhs.numberLinks_; 125 125 sosType_ = rhs.sosType_; 126 126 if (numberMembers_) { 127 weights_ = CoinCopyOfArray(rhs.weights_, numberMembers_);128 which_ = CoinCopyOfArray(rhs.which_, numberMembers_*numberLinks_);127 weights_ = CoinCopyOfArray(rhs.weights_, numberMembers_); 128 which_ = CoinCopyOfArray(rhs.which_, numberMembers_ * numberLinks_); 129 129 } else { 130 130 weights_ = NULL; … … 135 135 } 136 136 137 // Destructor 138 CbcLink::~CbcLink 139 { 140 delete 141 delete 137 // Destructor 138 CbcLink::~CbcLink() 139 { 140 delete[] weights_; 141 delete[] which_; 142 142 } 143 143 144 144 // Infeasibility  large is 0.5 145 double 146 CbcLink::infeasibility(int & 145 double 146 CbcLink::infeasibility(int &preferredWay) const 147 147 { 148 148 int j; 149 int firstNonZero =1;149 int firstNonZero = 1; 150 150 int lastNonZero = 1; 151 OsiSolverInterface * 152 const double * 151 OsiSolverInterface *solver = model_>solver(); 152 const double *solution = model_>testSolution(); 153 153 //const double * lower = solver>getColLower(); 154 const double * upper = solver>getColUpper(); 155 double integerTolerance = 156 model_>getDblParam(CbcModel::CbcIntegerTolerance); 154 const double *upper = solver>getColUpper(); 155 double integerTolerance = model_>getDblParam(CbcModel::CbcIntegerTolerance); 157 156 double weight = 0.0; 158 double sum = 0.0;157 double sum = 0.0; 159 158 160 159 // check bounds etc 161 double lastWeight =1.0e100;162 int base =0;163 for (j =0;j<numberMembers_;j++) {164 for (int k =0;k<numberLinks_;k++) {165 int iColumn = which_[base +k];160 double lastWeight = 1.0e100; 161 int base = 0; 162 for (j = 0; j < numberMembers_; j++) { 163 for (int k = 0; k < numberLinks_; k++) { 164 int iColumn = which_[base + k]; 166 165 //if (lower[iColumn]) 167 166 //throw CoinError("Non zero lower bound in CBCLink","infeasibility","CbcLink"); 168 if (lastWeight >=weights_[j]1.0e7)169 throw CoinError("Weights too close together in CBCLink", "infeasibility","CbcLink");170 double value = CoinMax(0.0, solution[iColumn]);167 if (lastWeight >= weights_[j]  1.0e7) 168 throw CoinError("Weights too close together in CBCLink", "infeasibility", "CbcLink"); 169 double value = CoinMax(0.0, solution[iColumn]); 171 170 sum += value; 172 if (value >integerTolerance&&upper[iColumn]) {171 if (value > integerTolerance && upper[iColumn]) { 173 172 // Possibly due to scaling a fixed variable might slip through 174 if (value >upper[iColumn]+1.0e8) {173 if (value > upper[iColumn] + 1.0e8) { 175 174 // Could change to #ifdef CBC_DEBUG 176 175 #ifndef NDEBUG 177 if (model_>messageHandler()>logLevel() >1)176 if (model_>messageHandler()>logLevel() > 1) 178 177 printf("** Variable %d (%d) has value %g and upper bound of %g\n", 179 iColumn,j,value,upper[iColumn]);178 iColumn, j, value, upper[iColumn]); 180 179 #endif 181 } 182 value = CoinMin(value,upper[iColumn]);183 weight += weights_[j] *value;184 if (firstNonZero <0)185 firstNonZero =j;186 lastNonZero =j;180 } 181 value = CoinMin(value, upper[iColumn]); 182 weight += weights_[j] * value; 183 if (firstNonZero < 0) 184 firstNonZero = j; 185 lastNonZero = j; 187 186 } 188 187 } … … 190 189 } 191 190 double valueInfeasibility; 192 preferredWay =1;193 if (lastNonZero firstNonZero>=sosType_) {191 preferredWay = 1; 192 if (lastNonZero  firstNonZero >= sosType_) { 194 193 // find where to branch 195 assert (sum>0.0);194 assert(sum > 0.0); 196 195 weight /= sum; 197 valueInfeasibility = lastNonZero firstNonZero+1;198 valueInfeasibility *= 0.5 /((double)numberMembers_);196 valueInfeasibility = lastNonZero  firstNonZero + 1; 197 valueInfeasibility *= 0.5 / ((double)numberMembers_); 199 198 //#define DISTANCE 200 199 #ifdef DISTANCE 201 assert (sosType_==1); // code up200 assert(sosType_ == 1); // code up 202 201 /* may still be satisfied. 203 202 For LOS type 2 we might wish to move coding around … … 205 204 */ 206 205 int iWhere; 207 bool possible =false;208 for (iWhere =firstNonZero;iWhere<=lastNonZero;iWhere++) {209 if (fabs(weight weights_[iWhere])<1.0e8) {210 possible=true;211 206 bool possible = false; 207 for (iWhere = firstNonZero; iWhere <= lastNonZero; iWhere++) { 208 if (fabs(weight  weights_[iWhere]) < 1.0e8) { 209 possible = true; 210 break; 212 211 } 213 212 } 214 213 if (possible) { 215 214 // One could move some of this (+ arrays) into model_ 216 const CoinPackedMatrix * 217 const double * 218 const int * 219 const CoinBigIndex * 220 const int * 221 const double * 222 const double * 223 const double * 215 const CoinPackedMatrix *matrix = solver>getMatrixByCol(); 216 const double *element = matrix>getMutableElements(); 217 const int *row = matrix>getIndices(); 218 const CoinBigIndex *columnStart = matrix>getVectorStarts(); 219 const int *columnLength = matrix>getVectorLengths(); 220 const double *rowSolution = solver>getRowActivity(); 221 const double *rowLower = solver>getRowLower(); 222 const double *rowUpper = solver>getRowUpper(); 224 223 int numberRows = matrix>getNumRows(); 225 double * array = new double[numberRows];226 CoinZeroN(array, numberRows);227 int * which = new int[numberRows];228 int n =0;229 int base =numberLinks_*firstNonZero;230 for (j =firstNonZero;j<=lastNonZero;j++) {231 for (int k=0;k<numberLinks_;k++) {232 int iColumn = which_[base+k];233 double value = CoinMax(0.0,solution[iColumn]);234 if (value>integerTolerance&&upper[iColumn]) {235 value = CoinMin(value,upper[iColumn]);236 for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {237 238 239 240 a += value*element[j];241 242 243 244 which[n++]=iRow;245 a=value*element[j];246 assert(a);247 248 array[iRow]=a;249 250 251 252 253 } 254 base =numberLinks_*iWhere;255 for (int k =0;k<numberLinks_;k++) {256 int iColumn = which_[base+k];257 258 for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {259 260 261 262 a = value*element[j];263 264 265 266 which[n++]=iRow;267 a=value*element[j];268 assert(a);269 270 array[iRow]=a;271 272 } 273 for (j =0;j<n;j++) {274 275 276 277 if (distance>1.0e8) {278 if (distance+rowSolution[iRow]>rowUpper[iRow]+1.0e8) {279 possible=false;280 281 282 } else if (distance<1.0e8) {283 if (distance+rowSolution[iRow]<rowLower[iRow]1.0e8) {284 possible=false;285 286 } 287 288 } 289 for (j =0;j<n;j++)290 array[which[j]]=0.0;291 delete 292 delete 224 double *array = new double[numberRows]; 225 CoinZeroN(array, numberRows); 226 int *which = new int[numberRows]; 227 int n = 0; 228 int base = numberLinks_ * firstNonZero; 229 for (j = firstNonZero; j <= lastNonZero; j++) { 230 for (int k = 0; k < numberLinks_; k++) { 231 int iColumn = which_[base + k]; 232 double value = CoinMax(0.0, solution[iColumn]); 233 if (value > integerTolerance && upper[iColumn]) { 234 value = CoinMin(value, upper[iColumn]); 235 for (int j = columnStart[iColumn]; j < columnStart[iColumn] + columnLength[iColumn]; j++) { 236 int iRow = row[j]; 237 double a = array[iRow]; 238 if (a) { 239 a += value * element[j]; 240 if (!a) 241 a = 1.0e100; 242 } else { 243 which[n++] = iRow; 244 a = value * element[j]; 245 assert(a); 246 } 247 array[iRow] = a; 248 } 249 } 250 } 251 base += numberLinks_; 252 } 253 base = numberLinks_ * iWhere; 254 for (int k = 0; k < numberLinks_; k++) { 255 int iColumn = which_[base + k]; 256 const double value = 1.0; 257 for (int j = columnStart[iColumn]; j < columnStart[iColumn] + columnLength[iColumn]; j++) { 258 int iRow = row[j]; 259 double a = array[iRow]; 260 if (a) { 261 a = value * element[j]; 262 if (!a) 263 a = 1.0e100; 264 } else { 265 which[n++] = iRow; 266 a = value * element[j]; 267 assert(a); 268 } 269 array[iRow] = a; 270 } 271 } 272 for (j = 0; j < n; j++) { 273 int iRow = which[j]; 274 // moving to point will increase row solution by this 275 double distance = array[iRow]; 276 if (distance > 1.0e8) { 277 if (distance + rowSolution[iRow] > rowUpper[iRow] + 1.0e8) { 278 possible = false; 279 break; 280 } 281 } else if (distance < 1.0e8) { 282 if (distance + rowSolution[iRow] < rowLower[iRow]  1.0e8) { 283 possible = false; 284 break; 285 } 286 } 287 } 288 for (j = 0; j < n; j++) 289 array[which[j]] = 0.0; 290 delete[] array; 291 delete[] which; 293 292 if (possible) { 294 valueInfeasibility=0.0;295 printf("possible %d %d %d\n",firstNonZero,lastNonZero,iWhere);293 valueInfeasibility = 0.0; 294 printf("possible %d %d %d\n", firstNonZero, lastNonZero, iWhere); 296 295 } 297 296 } … … 304 303 305 304 // This looks at solution and sets bounds to contain solution 306 void 307 CbcLink::feasibleRegion() 305 void CbcLink::feasibleRegion() 308 306 { 309 307 int j; 310 int firstNonZero =1;308 int firstNonZero = 1; 311 309 int lastNonZero = 1; 312 OsiSolverInterface * solver = model_>solver(); 313 const double * solution = model_>testSolution(); 314 const double * upper = solver>getColUpper(); 315 double integerTolerance = 316 model_>getDblParam(CbcModel::CbcIntegerTolerance); 310 OsiSolverInterface *solver = model_>solver(); 311 const double *solution = model_>testSolution(); 312 const double *upper = solver>getColUpper(); 313 double integerTolerance = model_>getDblParam(CbcModel::CbcIntegerTolerance); 317 314 double weight = 0.0; 318 double sum = 0.0;319 320 int base =0;321 for (j =0;j<numberMembers_;j++) {322 for (int k =0;k<numberLinks_;k++) {323 int iColumn = which_[base +k];324 double value = CoinMax(0.0, solution[iColumn]);315 double sum = 0.0; 316 317 int base = 0; 318 for (j = 0; j < numberMembers_; j++) { 319 for (int k = 0; k < numberLinks_; k++) { 320 int iColumn = which_[base + k]; 321 double value = CoinMax(0.0, solution[iColumn]); 325 322 sum += value; 326 if (value >integerTolerance&&upper[iColumn]) {327 weight += weights_[j] *value;328 if (firstNonZero <0)329 firstNonZero =j;330 lastNonZero =j;323 if (value > integerTolerance && upper[iColumn]) { 324 weight += weights_[j] * value; 325 if (firstNonZero < 0) 326 firstNonZero = j; 327 lastNonZero = j; 331 328 } 332 329 } … … 334 331 } 335 332 #ifdef DISTANCE 336 if (lastNonZero firstNonZero>sosType_1) {333 if (lastNonZero  firstNonZero > sosType_  1) { 337 334 /* may still be satisfied. 338 335 For LOS type 2 we might wish to move coding around … … 340 337 */ 341 338 int iWhere; 342 bool possible =false;343 for (iWhere =firstNonZero;iWhere<=lastNonZero;iWhere++) {344 if (fabs(weight weights_[iWhere])<1.0e8) {345 possible=true;346 339 bool possible = false; 340 for (iWhere = firstNonZero; iWhere <= lastNonZero; iWhere++) { 341 if (fabs(weight  weights_[iWhere]) < 1.0e8) { 342 possible = true; 343 break; 347 344 } 348 345 } 349 346 if (possible) { 350 347 // One could move some of this (+ arrays) into model_ 351 const CoinPackedMatrix * 352 const double * 353 const int * 354 const CoinBigIndex * 355 const int * 356 const double * 357 const double * 358 const double * 348 const CoinPackedMatrix *matrix = solver>getMatrixByCol(); 349 const double *element = matrix>getMutableElements(); 350 const int *row = matrix>getIndices(); 351 const CoinBigIndex *columnStart = matrix>getVectorStarts(); 352 const int *columnLength = matrix>getVectorLengths(); 353 const double *rowSolution = solver>getRowActivity(); 354 const double *rowLower = solver>getRowLower(); 355 const double *rowUpper = solver>getRowUpper(); 359 356 int numberRows = matrix>getNumRows(); 360 double * array = new double[numberRows];361 CoinZeroN(array, numberRows);362 int * which = new int[numberRows];363 int n =0;364 int base =numberLinks_*firstNonZero;365 for (j =firstNonZero;j<=lastNonZero;j++) {366 for (int k=0;k<numberLinks_;k++) {367 int iColumn = which_[base+k];368 double value = CoinMax(0.0,solution[iColumn]);369 if (value>integerTolerance&&upper[iColumn]) {370 value = CoinMin(value,upper[iColumn]);371 for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {372 373 374 375 a += value*element[j];376 377 378 379 which[n++]=iRow;380 a=value*element[j];381 assert(a);382 383 array[iRow]=a;384 385 386 387 388 } 389 base =numberLinks_*iWhere;390 for (int k =0;k<numberLinks_;k++) {391 int iColumn = which_[base+k];392 393 for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {394 395 396 397 a = value*element[j];398 399 400 401 which[n++]=iRow;402 a=value*element[j];403 assert(a);404 405 array[iRow]=a;406 407 } 408 for (j =0;j<n;j++) {409 410 411 412 if (distance>1.0e8) {413 if (distance+rowSolution[iRow]>rowUpper[iRow]+1.0e8) {414 possible=false;415 416 417 } else if (distance<1.0e8) {418 if (distance+rowSolution[iRow]<rowLower[iRow]1.0e8) {419 possible=false;420 421 } 422 423 } 424 for (j =0;j<n;j++)425 array[which[j]]=0.0;426 delete 427 delete 357 double *array = new double[numberRows]; 358 CoinZeroN(array, numberRows); 359 int *which = new int[numberRows]; 360 int n = 0; 361 int base = numberLinks_ * firstNonZero; 362 for (j = firstNonZero; j <= lastNonZero; j++) { 363 for (int k = 0; k < numberLinks_; k++) { 364 int iColumn = which_[base + k]; 365 double value = CoinMax(0.0, solution[iColumn]); 366 if (value > integerTolerance && upper[iColumn]) { 367 value = CoinMin(value, upper[iColumn]); 368 for (int j = columnStart[iColumn]; j < columnStart[iColumn] + columnLength[iColumn]; j++) { 369 int iRow = row[j]; 370 double a = array[iRow]; 371 if (a) { 372 a += value * element[j]; 373 if (!a) 374 a = 1.0e100; 375 } else { 376 which[n++] = iRow; 377 a = value * element[j]; 378 assert(a); 379 } 380 array[iRow] = a; 381 } 382 } 383 } 384 base += numberLinks_; 385 } 386 base = numberLinks_ * iWhere; 387 for (int k = 0; k < numberLinks_; k++) { 388 int iColumn = which_[base + k]; 389 const double value = 1.0; 390 for (int j = columnStart[iColumn]; j < columnStart[iColumn] + columnLength[iColumn]; j++) { 391 int iRow = row[j]; 392 double a = array[iRow]; 393 if (a) { 394 a = value * element[j]; 395 if (!a) 396 a = 1.0e100; 397 } else { 398 which[n++] = iRow; 399 a = value * element[j]; 400 assert(a); 401 } 402 array[iRow] = a; 403 } 404 } 405 for (j = 0; j < n; j++) { 406 int iRow = which[j]; 407 // moving to point will increase row solution by this 408 double distance = array[iRow]; 409 if (distance > 1.0e8) { 410 if (distance + rowSolution[iRow] > rowUpper[iRow] + 1.0e8) { 411 possible = false; 412 break; 413 } 414 } else if (distance < 1.0e8) { 415 if (distance + rowSolution[iRow] < rowLower[iRow]  1.0e8) { 416 possible = false; 417 break; 418 } 419 } 420 } 421 for (j = 0; j < n; j++) 422 array[which[j]] = 0.0; 423 delete[] array; 424 delete[] which; 428 425 if (possible) { 429 printf("possible feas region %d %d %d\n",firstNonZero,lastNonZero,iWhere);430 firstNonZero=iWhere;431 lastNonZero=iWhere;426 printf("possible feas region %d %d %d\n", firstNonZero, lastNonZero, iWhere); 427 firstNonZero = iWhere; 428 lastNonZero = iWhere; 432 429 } 433 430 } 434 431 } 435 432 #else 436 assert (lastNonZerofirstNonZero<sosType_);433 assert(lastNonZero  firstNonZero < sosType_); 437 434 #endif 438 base =0;439 for (j =0;j<firstNonZero;j++) {440 for (int k =0;k<numberLinks_;k++) {441 int iColumn = which_[base +k];442 solver>setColUpper(iColumn, 0.0);435 base = 0; 436 for (j = 0; j < firstNonZero; j++) { 437 for (int k = 0; k < numberLinks_; k++) { 438 int iColumn = which_[base + k]; 439 solver>setColUpper(iColumn, 0.0); 443 440 } 444 441 base += numberLinks_; … … 446 443 // skip 447 444 base += numberLinks_; 448 for (j =lastNonZero+1;j<numberMembers_;j++) {449 for (int k =0;k<numberLinks_;k++) {450 int iColumn = which_[base +k];451 solver>setColUpper(iColumn, 0.0);445 for (j = lastNonZero + 1; j < numberMembers_; j++) { 446 for (int k = 0; k < numberLinks_; k++) { 447 int iColumn = which_[base + k]; 448 solver>setColUpper(iColumn, 0.0); 452 449 } 453 450 base += numberLinks_; … … 455 452 } 456 453 457 458 454 // Creates a branching object 459 CbcBranchingObject * 460 CbcLink::createCbcBranch(OsiSolverInterface * /*solver*/, const OsiBranchingInformation * /*info*/, int way) 455 CbcBranchingObject * 456 CbcLink::createCbcBranch(OsiSolverInterface * /*solver*/, const OsiBranchingInformation * /*info*/, int way) 461 457 { 462 458 int j; 463 const double * solution = model_>testSolution(); 464 double integerTolerance = 465 model_>getDblParam(CbcModel::CbcIntegerTolerance); 466 OsiSolverInterface * solver = model_>solver(); 467 const double * upper = solver>getColUpper(); 468 int firstNonFixed=1; 469 int lastNonFixed=1; 470 int firstNonZero=1; 459 const double *solution = model_>testSolution(); 460 double integerTolerance = model_>getDblParam(CbcModel::CbcIntegerTolerance); 461 OsiSolverInterface *solver = model_>solver(); 462 const double *upper = solver>getColUpper(); 463 int firstNonFixed = 1; 464 int lastNonFixed = 1; 465 int firstNonZero = 1; 471 466 int lastNonZero = 1; 472 467 double weight = 0.0; 473 double sum = 0.0;474 int base =0;475 for (j =0;j<numberMembers_;j++) {476 for (int k =0;k<numberLinks_;k++) {477 int iColumn = which_[base +k];468 double sum = 0.0; 469 int base = 0; 470 for (j = 0; j < numberMembers_; j++) { 471 for (int k = 0; k < numberLinks_; k++) { 472 int iColumn = which_[base + k]; 478 473 if (upper[iColumn]) { 479 double value = CoinMax(0.0, solution[iColumn]);474 double value = CoinMax(0.0, solution[iColumn]); 480 475 sum += value; 481 if (firstNonFixed <0)482 firstNonFixed =j;483 lastNonFixed =j;484 if (value >integerTolerance) {485 weight += weights_[j] *value;486 if (firstNonZero <0)487 firstNonZero =j;488 lastNonZero =j;476 if (firstNonFixed < 0) 477 firstNonFixed = j; 478 lastNonFixed = j; 479 if (value > integerTolerance) { 480 weight += weights_[j] * value; 481 if (firstNonZero < 0) 482 firstNonZero = j; 483 lastNonZero = j; 489 484 } 490 485 } … … 492 487 base += numberLinks_; 493 488 } 494 assert (lastNonZerofirstNonZero>=sosType_);489 assert(lastNonZero  firstNonZero >= sosType_); 495 490 // find where to branch 496 assert (sum>0.0);491 assert(sum > 0.0); 497 492 weight /= sum; 498 493 int iWhere; 499 double separator =0.0;500 for (iWhere =firstNonZero;iWhere<lastNonZero;iWhere++)501 if (weight <weights_[iWhere+1])494 double separator = 0.0; 495 for (iWhere = firstNonZero; iWhere < lastNonZero; iWhere++) 496 if (weight < weights_[iWhere + 1]) 502 497 break; 503 if (sosType_ ==1) {498 if (sosType_ == 1) { 504 499 // SOS 1 505 separator = 0.5 * (weights_[iWhere]+weights_[iWhere+1]);500 separator = 0.5 * (weights_[iWhere] + weights_[iWhere + 1]); 506 501 } else { 507 502 // SOS 2 508 if (iWhere==firstNonFixed) 509 iWhere++;; 510 if (iWhere==lastNonFixed1) 511 iWhere = lastNonFixed2; 512 separator = weights_[iWhere+1]; 503 if (iWhere == firstNonFixed) 504 iWhere++; 505 ; 506 if (iWhere == lastNonFixed  1) 507 iWhere = lastNonFixed  2; 508 separator = weights_[iWhere + 1]; 513 509 } 514 510 // create object 515 CbcBranchingObject * 516 branch = new CbcLinkBranchingObject(model_, this,way,separator);511 CbcBranchingObject *branch; 512 branch = new CbcLinkBranchingObject(model_, this, way, separator); 517 513 return branch; 518 514 } 519 515 // Useful constructor 520 CbcLinkBranchingObject::CbcLinkBranchingObject (CbcModel *model,521 const CbcLink *set,522 int way,523 524 : CbcBranchingObject(model,set>id(),way,0.5)516 CbcLinkBranchingObject::CbcLinkBranchingObject(CbcModel *model, 517 const CbcLink *set, 518 int way, 519 double separator) 520 : CbcBranchingObject(model, set>id(), way, 0.5) 525 521 { 526 522 set_ = set; … … 528 524 } 529 525 530 // Copy constructor 531 CbcLinkBranchingObject::CbcLinkBranchingObject ( const CbcLinkBranchingObject & rhs) :CbcBranchingObject(rhs) 532 { 533 set_=rhs.set_; 526 // Copy constructor 527 CbcLinkBranchingObject::CbcLinkBranchingObject(const CbcLinkBranchingObject &rhs) 528 : CbcBranchingObject(rhs) 529 { 530 set_ = rhs.set_; 534 531 separator_ = rhs.separator_; 535 532 } 536 533 537 // Assignment operator 538 CbcLinkBranchingObject & 539 CbcLinkBranchingObject::operator=( const CbcLinkBranchingObject&rhs)534 // Assignment operator 535 CbcLinkBranchingObject & 536 CbcLinkBranchingObject::operator=(const CbcLinkBranchingObject &rhs) 540 537 { 541 538 if (this != &rhs) { 542 539 CbcBranchingObject::operator=(rhs); 543 set_ =rhs.set_;540 set_ = rhs.set_; 544 541 separator_ = rhs.separator_; 545 542 } 546 543 return *this; 547 544 } 548 CbcBranchingObject * 545 CbcBranchingObject * 549 546 CbcLinkBranchingObject::clone() const 550 { 547 { 551 548 return (new CbcLinkBranchingObject(*this)); 552 549 } 553 550 554 555 // Destructor 556 CbcLinkBranchingObject::~CbcLinkBranchingObject () 551 // Destructor 552 CbcLinkBranchingObject::~CbcLinkBranchingObject() 557 553 { 558 554 } … … 563 559 int numberMembers = set_>numberMembers(); 564 560 int numberLinks = set_>numberLinks(); 565 const double * 566 const int * 567 OsiSolverInterface * 561 const double *weights = set_>weights(); 562 const int *which = set_>which(); 563 OsiSolverInterface *solver = model_>solver(); 568 564 //const double * lower = solver>getColLower(); 569 565 //const double * upper = solver>getColUpper(); 570 566 // *** for way  up means fix all those in down section 571 if (way_ <0) {567 if (way_ < 0) { 572 568 int i; 573 for ( i=0;i<numberMembers;i++) {569 for (i = 0; i < numberMembers; i++) { 574 570 if (weights[i] > separator_) 575 break; 576 } 577 assert (i<numberMembers); 578 int base=i*numberLinks;; 579 for (;i<numberMembers;i++) { 580 for (int k=0;k<numberLinks;k++) { 581 int iColumn = which[base+k]; 582 solver>setColUpper(iColumn,0.0); 571 break; 572 } 573 assert(i < numberMembers); 574 int base = i * numberLinks; 575 ; 576 for (; i < numberMembers; i++) { 577 for (int k = 0; k < numberLinks; k++) { 578 int iColumn = which[base + k]; 579 solver>setColUpper(iColumn, 0.0); 583 580 } 584 581 base += numberLinks; 585 582 } 586 way_ =1;// Swap direction583 way_ = 1; // Swap direction 587 584 } else { 588 585 int i; 589 int base =0;590 for ( i=0;i<numberMembers;i++) {586 int base = 0; 587 for (i = 0; i < numberMembers; i++) { 591 588 if (weights[i] >= separator_) { 592 589 break; 593 590 } else { 594 for (int k =0;k<numberLinks;k++) {595 int iColumn = which[base+k];596 solver>setColUpper(iColumn, 0.0);591 for (int k = 0; k < numberLinks; k++) { 592 int iColumn = which[base + k]; 593 solver>setColUpper(iColumn, 0.0); 597 594 } 598 595 base += numberLinks; 599 596 } 600 597 } 601 assert (i<numberMembers);602 way_ =1;// Swap direction598 assert(i < numberMembers); 599 way_ = 1; // Swap direction 603 600 } 604 601 return 0.0; 605 602 } 606 // Print what would happen 607 void 608 CbcLinkBranchingObject::print() 603 // Print what would happen 604 void CbcLinkBranchingObject::print() 609 605 { 610 606 int numberMembers = set_>numberMembers(); 611 607 int numberLinks = set_>numberLinks(); 612 const double * 613 const int * 614 OsiSolverInterface * 615 const double * 616 int first =numberMembers;617 int last =1;618 int numberFixed =0;619 int numberOther =0;608 const double *weights = set_>weights(); 609 const int *which = set_>which(); 610 OsiSolverInterface *solver = model_>solver(); 611 const double *upper = solver>getColUpper(); 612 int first = numberMembers; 613 int last = 1; 614 int numberFixed = 0; 615 int numberOther = 0; 620 616 int i; 621 int base =0;622 for ( i=0;i<numberMembers;i++) {623 for (int k =0;k<numberLinks;k++) {624 int iColumn = which[base +k];617 int base = 0; 618 for (i = 0; i < numberMembers; i++) { 619 for (int k = 0; k < numberLinks; k++) { 620 int iColumn = which[base + k]; 625 621 double bound = upper[iColumn]; 626 622 if (bound) { 627 first = CoinMin(first, i);628 last = CoinMax(last, i);623 first = CoinMin(first, i); 624 last = CoinMax(last, i); 629 625 } 630 626 } … … 632 628 } 633 629 // *** for way  up means fix all those in down section 634 base =0;635 if (way_ <0) {630 base = 0; 631 if (way_ < 0) { 636 632 printf("SOS Down"); 637 for ( i=0;i<numberMembers;i++) {638 if (weights[i] > separator_) 639 640 for (int k =0;k<numberLinks;k++) {641 int iColumn = which[base +k];633 for (i = 0; i < numberMembers; i++) { 634 if (weights[i] > separator_) 635 break; 636 for (int k = 0; k < numberLinks; k++) { 637 int iColumn = which[base + k]; 642 638 double bound = upper[iColumn]; 643 639 if (bound) 644 640 numberOther++; 645 641 } 646 642 base += numberLinks; 647 643 } 648 assert (i<numberMembers);649 for (; i<numberMembers;i++) {650 for (int k =0;k<numberLinks;k++) {651 int iColumn = which[base +k];644 assert(i < numberMembers); 645 for (; i < numberMembers; i++) { 646 for (int k = 0; k < numberLinks; k++) { 647 int iColumn = which[base + k]; 652 648 double bound = upper[iColumn]; 653 649 if (bound) … … 658 654 } else { 659 655 printf("SOS Up"); 660 for ( i=0;i<numberMembers;i++) {656 for (i = 0; i < numberMembers; i++) { 661 657 if (weights[i] >= separator_) 662 663 for (int k =0;k<numberLinks;k++) {664 int iColumn = which[base +k];658 break; 659 for (int k = 0; k < numberLinks; k++) { 660 int iColumn = which[base + k]; 665 661 double bound = upper[iColumn]; 666 662 if (bound) 667 663 numberFixed++; 668 664 } 669 665 base += numberLinks; 670 666 } 671 assert (i<numberMembers);672 for (; i<numberMembers;i++) {673 for (int k =0;k<numberLinks;k++) {674 int iColumn = which[base +k];667 assert(i < numberMembers); 668 for (; i < numberMembers; i++) { 669 for (int k = 0; k < numberLinks; k++) { 670 int iColumn = which[base + k]; 675 671 double bound = upper[iColumn]; 676 672 if (bound) … … 680 676 } 681 677 } 682 assert ((numberFixed%numberLinks)==0);683 assert ((numberOther%numberLinks)==0);678 assert((numberFixed % numberLinks) == 0); 679 assert((numberOther % numberLinks) == 0); 684 680 printf("  at %g, free range %d (%g) => %d (%g), %d would be fixed, %d other way\n", 685 separator_,first,weights[first],last,weights[last],numberFixed/numberLinks,686 numberOther/numberLinks);681 separator_, first, weights[first], last, weights[last], numberFixed / numberLinks, 682 numberOther / numberLinks); 687 683 } 688 684 /** Compare the \c this with \c brObj. \c this and \c brObj must be os the … … 695 691 */ 696 692 CbcRangeCompare 697 CbcLinkBranchingObject::compareBranchingObject 698 (const CbcBranchingObject* brObj, const bool replaceIfOverlap) 693 CbcLinkBranchingObject::compareBranchingObject(const CbcBranchingObject *brObj, const bool replaceIfOverlap) 699 694 { 700 695 throw("must implement");
Note: See TracChangeset
for help on using the changeset viewer.