Changeset 2385 for trunk/Clp/src/ClpDynamicExampleMatrix.cpp
 Timestamp:
 Jan 6, 2019 2:43:06 PM (4 months ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/Clp/src/ClpDynamicExampleMatrix.cpp
r2271 r2385 26 26 // Default Constructor 27 27 // 28 ClpDynamicExampleMatrix::ClpDynamicExampleMatrix 29 : ClpDynamicMatrix(),30 numberColumns_(0),31 startColumnGen_(NULL),32 rowGen_(NULL),33 elementGen_(NULL),34 costGen_(NULL),35 fullStartGen_(NULL),36 dynamicStatusGen_(NULL),37 idGen_(NULL),38 columnLowerGen_(NULL),39 40 { 41 28 ClpDynamicExampleMatrix::ClpDynamicExampleMatrix() 29 : ClpDynamicMatrix() 30 , numberColumns_(0) 31 , startColumnGen_(NULL) 32 , rowGen_(NULL) 33 , elementGen_(NULL) 34 , costGen_(NULL) 35 , fullStartGen_(NULL) 36 , dynamicStatusGen_(NULL) 37 , idGen_(NULL) 38 , columnLowerGen_(NULL) 39 , columnUpperGen_(NULL) 40 { 41 setType(25); 42 42 } 43 43 … … 45 45 // Copy constructor 46 46 // 47 ClpDynamicExampleMatrix::ClpDynamicExampleMatrix (const ClpDynamicExampleMatrix & rhs) 48 : ClpDynamicMatrix(rhs) 49 { 50 numberColumns_ = rhs.numberColumns_; 51 startColumnGen_ = ClpCopyOfArray(rhs.startColumnGen_, numberColumns_ + 1); 52 CoinBigIndex numberElements = startColumnGen_[numberColumns_]; 53 rowGen_ = ClpCopyOfArray(rhs.rowGen_, numberElements);; 54 elementGen_ = ClpCopyOfArray(rhs.elementGen_, numberElements);; 55 costGen_ = ClpCopyOfArray(rhs.costGen_, numberColumns_); 56 fullStartGen_ = ClpCopyOfArray(rhs.fullStartGen_, numberSets_ + 1); 57 dynamicStatusGen_ = ClpCopyOfArray(rhs.dynamicStatusGen_, numberColumns_); 58 idGen_ = ClpCopyOfArray(rhs.idGen_, maximumGubColumns_); 59 columnLowerGen_ = ClpCopyOfArray(rhs.columnLowerGen_, numberColumns_); 60 columnUpperGen_ = ClpCopyOfArray(rhs.columnUpperGen_, numberColumns_); 47 ClpDynamicExampleMatrix::ClpDynamicExampleMatrix(const ClpDynamicExampleMatrix &rhs) 48 : ClpDynamicMatrix(rhs) 49 { 50 numberColumns_ = rhs.numberColumns_; 51 startColumnGen_ = ClpCopyOfArray(rhs.startColumnGen_, numberColumns_ + 1); 52 CoinBigIndex numberElements = startColumnGen_[numberColumns_]; 53 rowGen_ = ClpCopyOfArray(rhs.rowGen_, numberElements); 54 ; 55 elementGen_ = ClpCopyOfArray(rhs.elementGen_, numberElements); 56 ; 57 costGen_ = ClpCopyOfArray(rhs.costGen_, numberColumns_); 58 fullStartGen_ = ClpCopyOfArray(rhs.fullStartGen_, numberSets_ + 1); 59 dynamicStatusGen_ = ClpCopyOfArray(rhs.dynamicStatusGen_, numberColumns_); 60 idGen_ = ClpCopyOfArray(rhs.idGen_, maximumGubColumns_); 61 columnLowerGen_ = ClpCopyOfArray(rhs.columnLowerGen_, numberColumns_); 62 columnUpperGen_ = ClpCopyOfArray(rhs.columnUpperGen_, numberColumns_); 61 63 } 62 64 63 65 /* This is the real constructor*/ 64 ClpDynamicExampleMatrix::ClpDynamicExampleMatrix(ClpSimplex * 65 int numberGubColumns, const CoinBigIndex *starts,66 const double * lower, const double *upper,67 const CoinBigIndex * startColumn, const int *row,68 const double * element, const double *cost,69 const double * columnLower, const double *columnUpper,70 const unsigned char *status,71 const unsigned char *dynamicStatus,72 73 74 75 { 76 77 78 79 80 81 82 delete[] startSet_;83 startSet_ = new int[numberSets_];84 delete[] next_;85 next_ = new int[maximumGubColumns_];86 delete[] row_;87 delete[] element_;88 delete[] startColumn_;89 delete[] cost_;90 delete[] columnLower_;91 delete[] columnUpper_;92 delete[] dynamicStatus_;93 delete[] status_;94 delete[] id_;95 96 row_ = new int[maximumElements_];97 element_ = new double[maximumElements_];98 startColumn_ = new CoinBigIndex [maximumGubColumns_+1];99 100 101 102 103 dynamicStatus_ = new unsigned char [2*maximumGubColumns_];104 105 106 107 108 109 110 111 112 113 114 115 idGen_ = new int[maximumGubColumns_];116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 CoinSort_2(rowGen_ + startColumnGen_[i], rowGen_ + startColumnGen_[i+1], elementGen_ + startColumnGen_[i]);133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 double *columnUpper = model>columnUpper();157 for(i = firstDynamic_; i < lastDynamic_; i++)158 159 160 status_ = new unsigned char [2*numberSets_+4];161 162 memcpy(status_,status, numberSets_ * sizeof(char));163 assert(dynamicStatus);164 165 assert(numberIds);166 167 assert(!numberIds);168 169 170 171 172 173 174 dynamicStatusGen_ = new unsigned char[numberColumns_];175 176 177 178 179 180 181 182 183 184 addColumn(startColumnGen_[sequence+1]  start,185 186 187 188 189 190 191 192 193 194 195 196 int *set = new int[numberColumns_];197 198 for (CoinBigIndex j = fullStartGen_[iSet]; j < fullStartGen_[iSet+1]; j++)199 200 201 202 203 204 addColumn(startColumnGen_[sequence+1]  start,205 206 207 208 209 210 211 212 213 214 delete[] set;215 216 217 218 219 220 66 ClpDynamicExampleMatrix::ClpDynamicExampleMatrix(ClpSimplex *model, int numberSets, 67 int numberGubColumns, const CoinBigIndex *starts, 68 const double *lower, const double *upper, 69 const CoinBigIndex *startColumn, const int *row, 70 const double *element, const double *cost, 71 const double *columnLower, const double *columnUpper, 72 const unsigned char *status, 73 const unsigned char *dynamicStatus, 74 int numberIds, const int *ids) 75 : ClpDynamicMatrix(model, numberSets, 0, NULL, lower, upper, NULL, NULL, NULL, NULL, NULL, NULL, 76 NULL, NULL) 77 { 78 setType(25); 79 numberColumns_ = numberGubColumns; 80 // start with safe values  then experiment 81 maximumGubColumns_ = numberColumns_; 82 maximumElements_ = startColumn[numberColumns_]; 83 // delete odd stuff created by ClpDynamicMatrix constructor 84 delete[] startSet_; 85 startSet_ = new int[numberSets_]; 86 delete[] next_; 87 next_ = new int[maximumGubColumns_]; 88 delete[] row_; 89 delete[] element_; 90 delete[] startColumn_; 91 delete[] cost_; 92 delete[] columnLower_; 93 delete[] columnUpper_; 94 delete[] dynamicStatus_; 95 delete[] status_; 96 delete[] id_; 97 // and size correctly 98 row_ = new int[maximumElements_]; 99 element_ = new double[maximumElements_]; 100 startColumn_ = new CoinBigIndex[maximumGubColumns_ + 1]; 101 // say no columns yet 102 numberGubColumns_ = 0; 103 startColumn_[0] = 0; 104 cost_ = new double[maximumGubColumns_]; 105 dynamicStatus_ = new unsigned char[2 * maximumGubColumns_]; 106 memset(dynamicStatus_, 0, maximumGubColumns_); 107 id_ = new int[maximumGubColumns_]; 108 if (columnLower) 109 columnLower_ = new double[maximumGubColumns_]; 110 else 111 columnLower_ = NULL; 112 if (columnUpper) 113 columnUpper_ = new double[maximumGubColumns_]; 114 else 115 columnUpper_ = NULL; 116 // space for ids 117 idGen_ = new int[maximumGubColumns_]; 118 int iSet; 119 for (iSet = 0; iSet < numberSets_; iSet++) 120 startSet_[iSet] = 1; 121 // This starts code specific to this storage method 122 CoinBigIndex i; 123 fullStartGen_ = ClpCopyOfArray(starts, numberSets_ + 1); 124 startColumnGen_ = ClpCopyOfArray(startColumn, numberColumns_ + 1); 125 CoinBigIndex numberElements = startColumnGen_[numberColumns_]; 126 rowGen_ = ClpCopyOfArray(row, numberElements); 127 elementGen_ = new double[numberElements]; 128 for (i = 0; i < numberElements; i++) 129 elementGen_[i] = element[i]; 130 costGen_ = new double[numberColumns_]; 131 for (i = 0; i < numberColumns_; i++) { 132 costGen_[i] = cost[i]; 133 // I don't think I need sorted but ... 134 CoinSort_2(rowGen_ + startColumnGen_[i], rowGen_ + startColumnGen_[i + 1], elementGen_ + startColumnGen_[i]); 135 } 136 if (columnLower) { 137 columnLowerGen_ = new double[numberColumns_]; 138 for (i = 0; i < numberColumns_; i++) { 139 columnLowerGen_[i] = columnLower[i]; 140 if (columnLowerGen_[i]) { 141 printf("Nonzero lower bounds not allowed  subtract from model\n"); 142 abort(); 143 } 144 } 145 } else { 146 columnLowerGen_ = NULL; 147 } 148 if (columnUpper) { 149 columnUpperGen_ = new double[numberColumns_]; 150 for (i = 0; i < numberColumns_; i++) 151 columnUpperGen_[i] = columnUpper[i]; 152 } else { 153 columnUpperGen_ = NULL; 154 } 155 // end specific coding 156 if (columnUpper_) { 157 // set all upper bounds so we have enough space 158 double *columnUpper = model>columnUpper(); 159 for (i = firstDynamic_; i < lastDynamic_; i++) 160 columnUpper[i] = 1.0e10; 161 } 162 status_ = new unsigned char[2 * numberSets_ + 4]; 163 if (status) { 164 memcpy(status_, status, numberSets_ * sizeof(char)); 165 assert(dynamicStatus); 166 CoinMemcpyN(dynamicStatus, numberIds, dynamicStatus_); 167 assert(numberIds); 168 } else { 169 assert(!numberIds); 170 memset(status_, 0, numberSets_); 171 for (int i = 0; i < numberSets_; i++) { 172 // make slack key 173 setStatus(i, ClpSimplex::basic); 174 } 175 } 176 dynamicStatusGen_ = new unsigned char[numberColumns_]; 177 memset(dynamicStatusGen_, 0, numberColumns_); // for clarity 178 for (int i = 0; i < numberColumns_; i++) 179 setDynamicStatusGen(i, atLowerBound); 180 // Populate with enough columns 181 if (!numberIds) { 182 // This could be made more sophisticated 183 for (iSet = 0; iSet < numberSets_; iSet++) { 184 CoinBigIndex sequence = fullStartGen_[iSet]; 185 CoinBigIndex start = startColumnGen_[sequence]; 186 addColumn(startColumnGen_[sequence + 1]  start, 187 rowGen_ + start, 188 elementGen_ + start, 189 costGen_[sequence], 190 columnLowerGen_ ? columnLowerGen_[sequence] : 0, 191 columnUpperGen_ ? columnUpperGen_[sequence] : 1.0e30, 192 iSet, getDynamicStatusGen(sequence)); 193 idGen_[iSet] = sequence; // say which one in 194 setDynamicStatusGen(sequence, inSmall); 195 } 196 } else { 197 // put back old ones 198 int *set = new int[numberColumns_]; 199 for (iSet = 0; iSet < numberSets_; iSet++) { 200 for (CoinBigIndex j = fullStartGen_[iSet]; j < fullStartGen_[iSet + 1]; j++) 201 set[j] = iSet; 202 } 203 for (int i = 0; i < numberIds; i++) { 204 int sequence = ids[i]; 205 CoinBigIndex start = startColumnGen_[sequence]; 206 addColumn(startColumnGen_[sequence + 1]  start, 207 rowGen_ + start, 208 elementGen_ + start, 209 costGen_[sequence], 210 columnLowerGen_ ? columnLowerGen_[sequence] : 0, 211 columnUpperGen_ ? columnUpperGen_[sequence] : 1.0e30, 212 set[sequence], getDynamicStatus(i)); 213 idGen_[iSet] = sequence; // say which one in 214 setDynamicStatusGen(sequence, inSmall); 215 } 216 delete[] set; 217 } 218 if (!status) { 219 gubCrash(); 220 } else { 221 initialProblem(); 222 } 221 223 } 222 224 #if 0 … … 378 380 // Destructor 379 381 // 380 ClpDynamicExampleMatrix::~ClpDynamicExampleMatrix 381 { 382 delete[] startColumnGen_;383 delete[] rowGen_;384 delete[] elementGen_;385 delete[] costGen_;386 delete[] fullStartGen_;387 delete[] dynamicStatusGen_;388 delete[] idGen_;389 delete[] columnLowerGen_;390 delete[] columnUpperGen_;382 ClpDynamicExampleMatrix::~ClpDynamicExampleMatrix() 383 { 384 delete[] startColumnGen_; 385 delete[] rowGen_; 386 delete[] elementGen_; 387 delete[] costGen_; 388 delete[] fullStartGen_; 389 delete[] dynamicStatusGen_; 390 delete[] idGen_; 391 delete[] columnLowerGen_; 392 delete[] columnUpperGen_; 391 393 } 392 394 … … 395 397 // 396 398 ClpDynamicExampleMatrix & 397 ClpDynamicExampleMatrix::operator=(const ClpDynamicExampleMatrix &rhs)398 { 399 400 401 402 delete[] startColumnGen_;403 delete[] rowGen_;404 delete[] elementGen_;405 delete[] costGen_;406 delete[] fullStartGen_;407 delete[] dynamicStatusGen_;408 delete[] idGen_;409 delete[] columnLowerGen_;410 delete[] columnUpperGen_;411 412 413 414 415 416 417 418 419 420 421 422 399 ClpDynamicExampleMatrix::operator=(const ClpDynamicExampleMatrix &rhs) 400 { 401 if (this != &rhs) { 402 ClpDynamicMatrix::operator=(rhs); 403 numberColumns_ = rhs.numberColumns_; 404 delete[] startColumnGen_; 405 delete[] rowGen_; 406 delete[] elementGen_; 407 delete[] costGen_; 408 delete[] fullStartGen_; 409 delete[] dynamicStatusGen_; 410 delete[] idGen_; 411 delete[] columnLowerGen_; 412 delete[] columnUpperGen_; 413 startColumnGen_ = ClpCopyOfArray(rhs.startColumnGen_, numberColumns_ + 1); 414 CoinBigIndex numberElements = startColumnGen_[numberColumns_]; 415 rowGen_ = ClpCopyOfArray(rhs.rowGen_, numberElements); 416 elementGen_ = ClpCopyOfArray(rhs.elementGen_, numberElements); 417 costGen_ = ClpCopyOfArray(rhs.costGen_, numberColumns_); 418 fullStartGen_ = ClpCopyOfArray(rhs.fullStartGen_, numberSets_ + 1); 419 dynamicStatusGen_ = ClpCopyOfArray(rhs.dynamicStatusGen_, numberColumns_); 420 idGen_ = ClpCopyOfArray(rhs.idGen_, maximumGubColumns_); 421 columnLowerGen_ = ClpCopyOfArray(rhs.columnLowerGen_, numberColumns_); 422 columnUpperGen_ = ClpCopyOfArray(rhs.columnUpperGen_, numberColumns_); 423 } 424 return *this; 423 425 } 424 426 // 425 427 // Clone 426 428 // 427 ClpMatrixBase * 428 { 429 429 ClpMatrixBase *ClpDynamicExampleMatrix::clone() const 430 { 431 return new ClpDynamicExampleMatrix(*this); 430 432 } 431 433 // Partial pricing 432 void 433 ClpDynamicExampleMatrix::partialPricing(ClpSimplex * model, double startFraction, double endFraction, 434 int & bestSequence, int & numberWanted) 435 { 436 numberWanted = currentWanted_; 437 assert(!model>rowScale()); 438 if (!numberSets_) { 439 // no gub 440 ClpPackedMatrix::partialPricing(model, startFraction, endFraction, bestSequence, numberWanted); 441 } else { 442 // and do some proportion of full set 443 int startG2 = static_cast<int> (startFraction * numberSets_); 444 int endG2 = static_cast<int> (endFraction * numberSets_ + 0.1); 445 endG2 = CoinMin(endG2, numberSets_); 446 //printf("gub price  set start %d end %d\n", 447 // startG2,endG2); 448 double tolerance = model>currentDualTolerance(); 449 double * reducedCost = model>djRegion(); 450 const double * duals = model>dualRowSolution(); 451 double bestDj; 452 int numberRows = model>numberRows(); 453 int slackOffset = lastDynamic_ + numberRows; 454 int structuralOffset = slackOffset + numberSets_; 455 int structuralOffset2 = structuralOffset + maximumGubColumns_; 456 // If nothing found yet can go all the way to end 457 int endAll = endG2; 458 if (bestSequence < 0 && !startG2) 459 endAll = numberSets_; 460 if (bestSequence >= 0) { 461 if (bestSequence != savedBestSequence_) 462 bestDj = fabs(reducedCost[bestSequence]); // dj from slacks or permanent 463 else 464 bestDj = savedBestDj_; 465 } else { 466 bestDj = tolerance; 467 } 468 int saveSequence = bestSequence; 469 double djMod = 0.0; 470 double bestDjMod = 0.0; 471 //printf("iteration %d start %d end %d  wanted %d\n",model>numberIterations(), 472 // startG2,endG2,numberWanted); 473 int bestSet = 1; 474 int minSet = minimumObjectsScan_ < 0 ? 5 : minimumObjectsScan_; 475 int minNeg = minimumGoodReducedCosts_ < 0 ? 5 : minimumGoodReducedCosts_; 476 for (int iSet = startG2; iSet < endAll; iSet++) { 477 if (numberWanted + minNeg < originalWanted_ && iSet > startG2 + minSet) { 478 // give up 479 numberWanted = 0; 480 break; 481 } else if (iSet == endG2 && bestSequence >= 0) { 482 break; 483 } 484 int gubRow = toIndex_[iSet]; 485 if (gubRow >= 0) { 486 djMod = duals[gubRow+numberStaticRows_]; // have I got sign right? 487 } else { 488 int iBasic = keyVariable_[iSet]; 489 if (iBasic >= numberColumns_) { 490 djMod = 0.0; // set not in 491 } else { 492 // get dj without 493 djMod = 0.0; 494 for (CoinBigIndex j = startColumn_[iBasic]; 495 j < startColumn_[iBasic+1]; j++) { 496 int jRow = row_[j]; 497 djMod = duals[jRow] * element_[j]; 498 } 499 djMod += cost_[iBasic]; 500 // See if gub slack possible  dj is djMod 501 if (getStatus(iSet) == ClpSimplex::atLowerBound) { 502 double value = djMod; 503 if (value > tolerance) { 504 numberWanted; 505 if (value > bestDj) { 506 // check flagged variable and correct dj 507 if (!flagged(iSet)) { 508 bestDj = value; 509 bestSequence = slackOffset + iSet; 510 bestDjMod = djMod; 511 bestSet = iSet; 512 } else { 513 // just to make sure we don't exit before got something 514 numberWanted++; 515 abort(); 516 } 517 } 518 } 519 } else if (getStatus(iSet) == ClpSimplex::atUpperBound) { 520 double value = djMod; 521 if (value > tolerance) { 522 numberWanted; 523 if (value > bestDj) { 524 // check flagged variable and correct dj 525 if (!flagged(iSet)) { 526 bestDj = value; 527 bestSequence = slackOffset + iSet; 528 bestDjMod = djMod; 529 bestSet = iSet; 530 } else { 531 // just to make sure we don't exit before got something 532 numberWanted++; 533 abort(); 534 } 535 } 536 } 537 } 538 } 539 } 540 // do ones in small 541 CoinBigIndex iSequence = startSet_[iSet]; 542 while (iSequence >= 0) { 543 DynamicStatus status = getDynamicStatus(iSequence); 544 if (status == atLowerBound  status == atUpperBound) { 545 double value = cost_[iSequence]  djMod; 546 for (CoinBigIndex j = startColumn_[iSequence]; 547 j < startColumn_[iSequence+1]; j++) { 548 int jRow = row_[j]; 549 value = duals[jRow] * element_[j]; 550 } 551 // change sign if at lower bound 552 if (status == atLowerBound) 553 value = value; 554 if (value > tolerance) { 555 numberWanted; 556 if (value > bestDj) { 557 // check flagged variable and correct dj 558 if (!flagged(iSequence)) { 559 bestDj = value; 560 bestSequence = structuralOffset + iSequence; 561 bestDjMod = djMod; 562 bestSet = iSet; 563 } else { 564 // just to make sure we don't exit before got something 565 numberWanted++; 566 } 567 } 568 } 569 } 570 iSequence = next_[iSequence]; //onto next in set 571 } 572 // and now get best by column generation 573 // If no upper bounds we may not need status test 574 for (iSequence = fullStartGen_[iSet]; iSequence < fullStartGen_[iSet+1]; iSequence++) { 575 DynamicStatus status = getDynamicStatusGen(iSequence); 576 assert (status != atUpperBound && status != soloKey); 577 if (status == atLowerBound) { 578 double value = costGen_[iSequence]  djMod; 579 for (CoinBigIndex j = startColumnGen_[iSequence]; 580 j < startColumnGen_[iSequence+1]; j++) { 581 int jRow = rowGen_[j]; 582 value = duals[jRow] * elementGen_[j]; 583 } 584 // change sign as at lower bound 585 value = value; 586 if (value > tolerance) { 587 numberWanted; 588 if (value > bestDj) { 589 // check flagged variable and correct dj 590 if (!flaggedGen(iSequence)) { 591 bestDj = value; 592 bestSequence = structuralOffset2 + iSequence; 593 bestDjMod = djMod; 594 bestSet = iSet; 595 } else { 596 // just to make sure we don't exit before got something 597 numberWanted++; 598 } 599 } 600 } 601 } 602 } 603 if (numberWanted <= 0) { 604 numberWanted = 0; 605 break; 606 } 607 } 608 if (bestSequence != saveSequence) { 609 savedBestGubDual_ = bestDjMod; 610 savedBestDj_ = bestDj; 611 savedBestSequence_ = bestSequence; 612 savedBestSet_ = bestSet; 613 } 614 // Do packed part before gub 615 // always??? 616 // Resize so just do to gub 617 numberActiveColumns_ = firstDynamic_; 618 int saveMinNeg = minimumGoodReducedCosts_; 619 if (bestSequence >= 0) 620 minimumGoodReducedCosts_ = 2; 621 currentWanted_ = numberWanted; 622 ClpPackedMatrix::partialPricing(model, startFraction, endFraction, bestSequence, numberWanted); 623 numberActiveColumns_ = matrix_>getNumCols(); 624 minimumGoodReducedCosts_ = saveMinNeg; 625 // See if may be finished 626 if (!startG2 && bestSequence < 0) 627 infeasibilityWeight_ = model_>infeasibilityCost(); 628 else if (bestSequence >= 0) 629 infeasibilityWeight_ = 1.0; 630 currentWanted_ = numberWanted; 631 } 434 void ClpDynamicExampleMatrix::partialPricing(ClpSimplex *model, double startFraction, double endFraction, 435 int &bestSequence, int &numberWanted) 436 { 437 numberWanted = currentWanted_; 438 assert(!model>rowScale()); 439 if (!numberSets_) { 440 // no gub 441 ClpPackedMatrix::partialPricing(model, startFraction, endFraction, bestSequence, numberWanted); 442 } else { 443 // and do some proportion of full set 444 int startG2 = static_cast< int >(startFraction * numberSets_); 445 int endG2 = static_cast< int >(endFraction * numberSets_ + 0.1); 446 endG2 = CoinMin(endG2, numberSets_); 447 //printf("gub price  set start %d end %d\n", 448 // startG2,endG2); 449 double tolerance = model>currentDualTolerance(); 450 double *reducedCost = model>djRegion(); 451 const double *duals = model>dualRowSolution(); 452 double bestDj; 453 int numberRows = model>numberRows(); 454 int slackOffset = lastDynamic_ + numberRows; 455 int structuralOffset = slackOffset + numberSets_; 456 int structuralOffset2 = structuralOffset + maximumGubColumns_; 457 // If nothing found yet can go all the way to end 458 int endAll = endG2; 459 if (bestSequence < 0 && !startG2) 460 endAll = numberSets_; 461 if (bestSequence >= 0) { 462 if (bestSequence != savedBestSequence_) 463 bestDj = fabs(reducedCost[bestSequence]); // dj from slacks or permanent 464 else 465 bestDj = savedBestDj_; 466 } else { 467 bestDj = tolerance; 468 } 469 int saveSequence = bestSequence; 470 double djMod = 0.0; 471 double bestDjMod = 0.0; 472 //printf("iteration %d start %d end %d  wanted %d\n",model>numberIterations(), 473 // startG2,endG2,numberWanted); 474 int bestSet = 1; 475 int minSet = minimumObjectsScan_ < 0 ? 5 : minimumObjectsScan_; 476 int minNeg = minimumGoodReducedCosts_ < 0 ? 5 : minimumGoodReducedCosts_; 477 for (int iSet = startG2; iSet < endAll; iSet++) { 478 if (numberWanted + minNeg < originalWanted_ && iSet > startG2 + minSet) { 479 // give up 480 numberWanted = 0; 481 break; 482 } else if (iSet == endG2 && bestSequence >= 0) { 483 break; 484 } 485 int gubRow = toIndex_[iSet]; 486 if (gubRow >= 0) { 487 djMod = duals[gubRow + numberStaticRows_]; // have I got sign right? 488 } else { 489 int iBasic = keyVariable_[iSet]; 490 if (iBasic >= numberColumns_) { 491 djMod = 0.0; // set not in 492 } else { 493 // get dj without 494 djMod = 0.0; 495 for (CoinBigIndex j = startColumn_[iBasic]; 496 j < startColumn_[iBasic + 1]; j++) { 497 int jRow = row_[j]; 498 djMod = duals[jRow] * element_[j]; 499 } 500 djMod += cost_[iBasic]; 501 // See if gub slack possible  dj is djMod 502 if (getStatus(iSet) == ClpSimplex::atLowerBound) { 503 double value = djMod; 504 if (value > tolerance) { 505 numberWanted; 506 if (value > bestDj) { 507 // check flagged variable and correct dj 508 if (!flagged(iSet)) { 509 bestDj = value; 510 bestSequence = slackOffset + iSet; 511 bestDjMod = djMod; 512 bestSet = iSet; 513 } else { 514 // just to make sure we don't exit before got something 515 numberWanted++; 516 abort(); 517 } 518 } 519 } 520 } else if (getStatus(iSet) == ClpSimplex::atUpperBound) { 521 double value = djMod; 522 if (value > tolerance) { 523 numberWanted; 524 if (value > bestDj) { 525 // check flagged variable and correct dj 526 if (!flagged(iSet)) { 527 bestDj = value; 528 bestSequence = slackOffset + iSet; 529 bestDjMod = djMod; 530 bestSet = iSet; 531 } else { 532 // just to make sure we don't exit before got something 533 numberWanted++; 534 abort(); 535 } 536 } 537 } 538 } 539 } 540 } 541 // do ones in small 542 CoinBigIndex iSequence = startSet_[iSet]; 543 while (iSequence >= 0) { 544 DynamicStatus status = getDynamicStatus(iSequence); 545 if (status == atLowerBound  status == atUpperBound) { 546 double value = cost_[iSequence]  djMod; 547 for (CoinBigIndex j = startColumn_[iSequence]; 548 j < startColumn_[iSequence + 1]; j++) { 549 int jRow = row_[j]; 550 value = duals[jRow] * element_[j]; 551 } 552 // change sign if at lower bound 553 if (status == atLowerBound) 554 value = value; 555 if (value > tolerance) { 556 numberWanted; 557 if (value > bestDj) { 558 // check flagged variable and correct dj 559 if (!flagged(iSequence)) { 560 bestDj = value; 561 bestSequence = structuralOffset + iSequence; 562 bestDjMod = djMod; 563 bestSet = iSet; 564 } else { 565 // just to make sure we don't exit before got something 566 numberWanted++; 567 } 568 } 569 } 570 } 571 iSequence = next_[iSequence]; //onto next in set 572 } 573 // and now get best by column generation 574 // If no upper bounds we may not need status test 575 for (iSequence = fullStartGen_[iSet]; iSequence < fullStartGen_[iSet + 1]; iSequence++) { 576 DynamicStatus status = getDynamicStatusGen(iSequence); 577 assert(status != atUpperBound && status != soloKey); 578 if (status == atLowerBound) { 579 double value = costGen_[iSequence]  djMod; 580 for (CoinBigIndex j = startColumnGen_[iSequence]; 581 j < startColumnGen_[iSequence + 1]; j++) { 582 int jRow = rowGen_[j]; 583 value = duals[jRow] * elementGen_[j]; 584 } 585 // change sign as at lower bound 586 value = value; 587 if (value > tolerance) { 588 numberWanted; 589 if (value > bestDj) { 590 // check flagged variable and correct dj 591 if (!flaggedGen(iSequence)) { 592 bestDj = value; 593 bestSequence = structuralOffset2 + iSequence; 594 bestDjMod = djMod; 595 bestSet = iSet; 596 } else { 597 // just to make sure we don't exit before got something 598 numberWanted++; 599 } 600 } 601 } 602 } 603 } 604 if (numberWanted <= 0) { 605 numberWanted = 0; 606 break; 607 } 608 } 609 if (bestSequence != saveSequence) { 610 savedBestGubDual_ = bestDjMod; 611 savedBestDj_ = bestDj; 612 savedBestSequence_ = bestSequence; 613 savedBestSet_ = bestSet; 614 } 615 // Do packed part before gub 616 // always??? 617 // Resize so just do to gub 618 numberActiveColumns_ = firstDynamic_; 619 int saveMinNeg = minimumGoodReducedCosts_; 620 if (bestSequence >= 0) 621 minimumGoodReducedCosts_ = 2; 622 currentWanted_ = numberWanted; 623 ClpPackedMatrix::partialPricing(model, startFraction, endFraction, bestSequence, numberWanted); 624 numberActiveColumns_ = matrix_>getNumCols(); 625 minimumGoodReducedCosts_ = saveMinNeg; 626 // See if may be finished 627 if (!startG2 && bestSequence < 0) 628 infeasibilityWeight_ = model_>infeasibilityCost(); 629 else if (bestSequence >= 0) 630 infeasibilityWeight_ = 1.0; 631 currentWanted_ = numberWanted; 632 } 632 633 } 633 634 /* Creates a variable. This is called after partial pricing and may modify matrix. 634 635 May update bestSequence. 635 636 */ 636 void 637 ClpDynamicExampleMatrix::createVariable(ClpSimplex * model, int & bestSequence) 638 { 639 int numberRows = model>numberRows(); 640 int slackOffset = lastDynamic_ + numberRows; 641 int structuralOffset = slackOffset + numberSets_; 642 int bestSequence2 = savedBestSequence_  structuralOffset; 643 if (bestSequence2 >= 0) { 644 // See if needs new 645 if (bestSequence2 >= maximumGubColumns_) { 646 bestSequence2 = maximumGubColumns_; 647 int sequence = addColumn(startColumnGen_[bestSequence2+1]  startColumnGen_[bestSequence2], 648 rowGen_ + startColumnGen_[bestSequence2], 649 elementGen_ + startColumnGen_[bestSequence2], 650 costGen_[bestSequence2], 651 columnLowerGen_ ? columnLowerGen_[bestSequence2] : 0, 652 columnUpperGen_ ? columnUpperGen_[bestSequence2] : 1.0e30, 653 savedBestSet_, getDynamicStatusGen(bestSequence2)); 654 savedBestSequence_ = structuralOffset + sequence; 655 idGen_[sequence] = bestSequence2; 656 setDynamicStatusGen(bestSequence2, inSmall); 657 } 658 } 659 ClpDynamicMatrix::createVariable(model, bestSequence/*, bestSequence2*/); 660 // clear for next iteration 661 savedBestSequence_ = 1; 637 void ClpDynamicExampleMatrix::createVariable(ClpSimplex *model, int &bestSequence) 638 { 639 int numberRows = model>numberRows(); 640 int slackOffset = lastDynamic_ + numberRows; 641 int structuralOffset = slackOffset + numberSets_; 642 int bestSequence2 = savedBestSequence_  structuralOffset; 643 if (bestSequence2 >= 0) { 644 // See if needs new 645 if (bestSequence2 >= maximumGubColumns_) { 646 bestSequence2 = maximumGubColumns_; 647 int sequence = addColumn(startColumnGen_[bestSequence2 + 1]  startColumnGen_[bestSequence2], 648 rowGen_ + startColumnGen_[bestSequence2], 649 elementGen_ + startColumnGen_[bestSequence2], 650 costGen_[bestSequence2], 651 columnLowerGen_ ? columnLowerGen_[bestSequence2] : 0, 652 columnUpperGen_ ? columnUpperGen_[bestSequence2] : 1.0e30, 653 savedBestSet_, getDynamicStatusGen(bestSequence2)); 654 savedBestSequence_ = structuralOffset + sequence; 655 idGen_[sequence] = bestSequence2; 656 setDynamicStatusGen(bestSequence2, inSmall); 657 } 658 } 659 ClpDynamicMatrix::createVariable(model, bestSequence /*, bestSequence2*/); 660 // clear for next iteration 661 savedBestSequence_ = 1; 662 662 } 663 663 /* If addColumn forces compression then this allows descendant to know what to do. … … 665 665 Entries at upper bound (really nonzero) never go out (at present). 666 666 */ 667 void 668 ClpDynamicExampleMatrix::packDown(const int * in, int numberToPack) 669 { 670 int put = 0; 671 for (int i = 0; i < numberToPack; i++) { 672 int id = idGen_[i]; 673 if (in[i] >= 0) { 674 // stays in 675 assert (put == in[i]); // true for now 676 idGen_[put++] = id; 677 } else { 678 // out to lower bound 679 setDynamicStatusGen(id, atLowerBound); 680 } 681 } 682 assert (put == numberGubColumns_); 683 } 667 void ClpDynamicExampleMatrix::packDown(const int *in, int numberToPack) 668 { 669 int put = 0; 670 for (int i = 0; i < numberToPack; i++) { 671 int id = idGen_[i]; 672 if (in[i] >= 0) { 673 // stays in 674 assert(put == in[i]); // true for now 675 idGen_[put++] = id; 676 } else { 677 // out to lower bound 678 setDynamicStatusGen(id, atLowerBound); 679 } 680 } 681 assert(put == numberGubColumns_); 682 } 683 684 /* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2 685 */
Note: See TracChangeset
for help on using the changeset viewer.