[1854] | 1 | /* $Id: CbcHeuristicPivotAndFix.cpp 1854 2013-01-28 00:02:55Z forrest $ */ |
---|
[1055] | 2 | // Copyright (C) 2008, International Business Machines |
---|
| 3 | // Corporation and others. All Rights Reserved. |
---|
[1573] | 4 | // This code is licensed under the terms of the Eclipse Public License (EPL). |
---|
| 5 | |
---|
[1055] | 6 | #if defined(_MSC_VER) |
---|
| 7 | // Turn off compiler warning about long names |
---|
| 8 | # pragma warning(disable:4786) |
---|
| 9 | #endif |
---|
| 10 | #include <cassert> |
---|
| 11 | #include <cstdlib> |
---|
| 12 | #include <cmath> |
---|
| 13 | #include <cfloat> |
---|
[1065] | 14 | #include <vector> |
---|
[1055] | 15 | |
---|
| 16 | #include "OsiSolverInterface.hpp" |
---|
| 17 | #include "CbcModel.hpp" |
---|
| 18 | #include "CbcMessage.hpp" |
---|
| 19 | #include "CbcHeuristicPivotAndFix.hpp" |
---|
[1065] | 20 | #include "OsiClpSolverInterface.hpp" |
---|
| 21 | #include "CoinTime.hpp" |
---|
[1055] | 22 | |
---|
[1065] | 23 | //#define FORNOW |
---|
| 24 | |
---|
[1055] | 25 | // Default Constructor |
---|
[1286] | 26 | CbcHeuristicPivotAndFix::CbcHeuristicPivotAndFix() |
---|
| 27 | : CbcHeuristic() |
---|
[1055] | 28 | { |
---|
| 29 | } |
---|
| 30 | |
---|
| 31 | // Constructor with model - assumed before cuts |
---|
| 32 | |
---|
| 33 | CbcHeuristicPivotAndFix::CbcHeuristicPivotAndFix(CbcModel & model) |
---|
[1286] | 34 | : CbcHeuristic(model) |
---|
[1055] | 35 | { |
---|
| 36 | } |
---|
| 37 | |
---|
[1286] | 38 | // Destructor |
---|
[1055] | 39 | CbcHeuristicPivotAndFix::~CbcHeuristicPivotAndFix () |
---|
| 40 | { |
---|
| 41 | } |
---|
| 42 | |
---|
| 43 | // Clone |
---|
| 44 | CbcHeuristic * |
---|
| 45 | CbcHeuristicPivotAndFix::clone() const |
---|
| 46 | { |
---|
[1286] | 47 | return new CbcHeuristicPivotAndFix(*this); |
---|
[1055] | 48 | } |
---|
| 49 | // Create C++ lines to get to current state |
---|
[1286] | 50 | void |
---|
| 51 | CbcHeuristicPivotAndFix::generateCpp( FILE * fp) |
---|
[1055] | 52 | { |
---|
[1286] | 53 | CbcHeuristicPivotAndFix other; |
---|
| 54 | fprintf(fp, "0#include \"CbcHeuristicPivotAndFix.hpp\"\n"); |
---|
| 55 | fprintf(fp, "3 CbcHeuristicPivotAndFix heuristicPFX(*cbcModel);\n"); |
---|
| 56 | CbcHeuristic::generateCpp(fp, "heuristicPFX"); |
---|
| 57 | fprintf(fp, "3 cbcModel->addHeuristic(&heuristicPFX);\n"); |
---|
[1055] | 58 | } |
---|
| 59 | |
---|
[1286] | 60 | // Copy constructor |
---|
[1055] | 61 | CbcHeuristicPivotAndFix::CbcHeuristicPivotAndFix(const CbcHeuristicPivotAndFix & rhs) |
---|
[1286] | 62 | : |
---|
| 63 | CbcHeuristic(rhs) |
---|
[1055] | 64 | { |
---|
| 65 | } |
---|
| 66 | |
---|
[1286] | 67 | // Assignment operator |
---|
| 68 | CbcHeuristicPivotAndFix & |
---|
| 69 | CbcHeuristicPivotAndFix::operator=( const CbcHeuristicPivotAndFix & rhs) |
---|
[1055] | 70 | { |
---|
[1286] | 71 | if (this != &rhs) { |
---|
| 72 | CbcHeuristic::operator=(rhs); |
---|
| 73 | } |
---|
| 74 | return *this; |
---|
[1055] | 75 | } |
---|
| 76 | |
---|
| 77 | // Resets stuff if model changes |
---|
[1286] | 78 | void |
---|
[1271] | 79 | CbcHeuristicPivotAndFix::resetModel(CbcModel * /*model*/) |
---|
[1055] | 80 | { |
---|
[1286] | 81 | //CbcHeuristic::resetModel(model); |
---|
[1055] | 82 | } |
---|
| 83 | /* |
---|
| 84 | Comments needed |
---|
| 85 | Returns 1 if solution, 0 if not */ |
---|
| 86 | int |
---|
[1271] | 87 | CbcHeuristicPivotAndFix::solution(double & /*solutionValue*/, |
---|
[1286] | 88 | double * /*betterSolution*/) |
---|
[1055] | 89 | { |
---|
| 90 | |
---|
[1286] | 91 | numCouldRun_++; // Todo: Ask JJHF what this for. |
---|
| 92 | std::cout << "Entering Pivot-and-Fix Heuristic" << std::endl; |
---|
[1065] | 93 | |
---|
| 94 | #ifdef FORNOW |
---|
[1286] | 95 | std::cout << "Lucky you! You're in the Pivot-and-Fix Heuristic" << std::endl; |
---|
| 96 | // The struct should be moved to member data |
---|
[1065] | 97 | |
---|
| 98 | |
---|
| 99 | |
---|
[1286] | 100 | typedef struct { |
---|
| 101 | int numberSolutions; |
---|
| 102 | int maximumSolutions; |
---|
| 103 | int numberColumns; |
---|
| 104 | double ** solution; |
---|
| 105 | int * numberUnsatisfied; |
---|
| 106 | } clpSolution; |
---|
[1065] | 107 | |
---|
[1286] | 108 | double start = CoinCpuTime(); |
---|
[1065] | 109 | |
---|
[1286] | 110 | OsiClpSolverInterface * clpSolverOriginal |
---|
[1065] | 111 | = dynamic_cast<OsiClpSolverInterface *> (model_->solver()); |
---|
[1286] | 112 | assert (clpSolverOriginal); |
---|
| 113 | OsiClpSolverInterface *clpSolver(clpSolverOriginal); |
---|
[1065] | 114 | |
---|
[1286] | 115 | ClpSimplex * simplex = clpSolver->getModelPtr(); |
---|
[1065] | 116 | |
---|
[1286] | 117 | // Initialize the structure holding the solutions |
---|
| 118 | clpSolution solutions; |
---|
| 119 | // Set typeStruct field of ClpTrustedData struct to one. |
---|
| 120 | // This tells Clp it's "Mahdi!" |
---|
| 121 | ClpTrustedData trustedSolutions; |
---|
| 122 | trustedSolutions.typeStruct = 1; |
---|
| 123 | trustedSolutions.data = &solutions; |
---|
| 124 | solutions.numberSolutions = 0; |
---|
| 125 | solutions.maximumSolutions = 0; |
---|
| 126 | solutions.numberColumns = simplex->numberColumns(); |
---|
| 127 | solutions.solution = NULL; |
---|
| 128 | solutions.numberUnsatisfied = NULL; |
---|
| 129 | simplex->setTrustedUserPointer(&trustedSolutions); |
---|
[1065] | 130 | |
---|
[1286] | 131 | // Solve from all slack to get some points |
---|
| 132 | simplex->allSlackBasis(); |
---|
| 133 | simplex->primal(); |
---|
[1065] | 134 | |
---|
[1286] | 135 | // ------------------------------------------------- |
---|
| 136 | // Get the problem information |
---|
| 137 | // - get the number of cols and rows |
---|
| 138 | int numCols = clpSolver->getNumCols(); |
---|
| 139 | int numRows = clpSolver->getNumRows(); |
---|
| 140 | |
---|
| 141 | // - get the right hand side of the rows |
---|
| 142 | const double * rhs = clpSolver->getRightHandSide(); |
---|
| 143 | |
---|
| 144 | // - find the integer variables |
---|
| 145 | bool * varClassInt = new bool[numCols]; |
---|
| 146 | int numInt = 0; |
---|
| 147 | for (int i = 0; i < numCols; i++) { |
---|
| 148 | if (clpSolver->isContinuous(i)) |
---|
| 149 | varClassInt[i] = 0; |
---|
| 150 | else { |
---|
| 151 | varClassInt[i] = 1; |
---|
| 152 | numInt++; |
---|
| 153 | } |
---|
[1065] | 154 | } |
---|
[1286] | 155 | |
---|
| 156 | // -Get the rows sense |
---|
| 157 | const char * rowSense; |
---|
| 158 | rowSense = clpSolver->getRowSense(); |
---|
| 159 | |
---|
| 160 | // -Get the objective coefficients |
---|
| 161 | const double *objCoefficients = clpSolver->getObjCoefficients(); |
---|
| 162 | double *originalObjCoeff = new double [numCols]; |
---|
| 163 | for (int i = 0; i < numCols; i++) |
---|
| 164 | originalObjCoeff[i] = objCoefficients[i]; |
---|
| 165 | |
---|
| 166 | // -Get the matrix of the problem |
---|
| 167 | double ** matrix = new double * [numRows]; |
---|
| 168 | for (int i = 0; i < numRows; i++) { |
---|
| 169 | matrix[i] = new double[numCols]; |
---|
| 170 | for (int j = 0; j < numCols; j++) |
---|
| 171 | matrix[i][j] = 0; |
---|
[1065] | 172 | } |
---|
[1286] | 173 | const CoinPackedMatrix* matrixByRow = clpSolver->getMatrixByRow(); |
---|
| 174 | const double * matrixElements = matrixByRow->getElements(); |
---|
| 175 | const int * matrixIndices = matrixByRow->getIndices(); |
---|
| 176 | const int * matrixStarts = matrixByRow->getVectorStarts(); |
---|
| 177 | for (int j = 0; j < numRows; j++) { |
---|
| 178 | for (int i = matrixStarts[j]; i < matrixStarts[j+1]; i++) { |
---|
| 179 | matrix[j][matrixIndices[i]] = matrixElements[i]; |
---|
| 180 | } |
---|
[1065] | 181 | } |
---|
| 182 | |
---|
[1286] | 183 | // The newObj is the randomly perturbed constraint used to find new |
---|
| 184 | // corner points |
---|
| 185 | double * newObj = new double [numCols]; |
---|
| 186 | |
---|
| 187 | // Set the random seed |
---|
| 188 | srand ( time(NULL) + 1); |
---|
| 189 | int randNum; |
---|
| 190 | |
---|
| 191 | // We're going to add a new row to the LP formulation |
---|
| 192 | // after finding each new solution. |
---|
| 193 | // Adding a new row requires the new elements and the new indices. |
---|
| 194 | // The elements are original objective function coefficients. |
---|
| 195 | // The indicies are the (dense) columns indices stored in addRowIndex. |
---|
| 196 | // The rhs is the value of the new solution stored in solutionValue. |
---|
| 197 | int * addRowIndex = new int[numCols]; |
---|
| 198 | for (int i = 0; i < numCols; i++) |
---|
| 199 | addRowIndex[i] = i; |
---|
| 200 | |
---|
| 201 | // The number of feasible solutions found by the PF heuristic. |
---|
| 202 | // This controls the return code of the solution() method. |
---|
| 203 | int numFeasibles = 0; |
---|
| 204 | |
---|
| 205 | // Shuffle the rows |
---|
| 206 | int * index = new int [numRows]; |
---|
| 207 | for (int i = 0; i < numRows; i++) |
---|
| 208 | index[i] = i; |
---|
| 209 | for (int i = 0; i < numRows; i++) { |
---|
| 210 | int temp = index[i]; |
---|
| 211 | int randNumTemp = i + (rand() % (numRows - i)); |
---|
| 212 | index[i] = index[randNumTemp]; |
---|
| 213 | index[randNumTemp] = temp; |
---|
[1065] | 214 | } |
---|
| 215 | |
---|
[1286] | 216 | // In the clpSolution struct, we store a lot of column solutions. |
---|
| 217 | // For each perturb objective, we store the solution from each |
---|
| 218 | // iteration of the LP solve. |
---|
| 219 | // For each perturb objective, we look at the collection of |
---|
| 220 | // solutions to do something extremly intelligent :-) |
---|
| 221 | // We could (and should..and will :-) wipe out the block of |
---|
| 222 | // solutions when we're done with them. But for now, we just move on |
---|
| 223 | // and store the next block of solutions for the next (perturbed) |
---|
| 224 | // objective. |
---|
| 225 | // The variable startIndex tells us where the new block begins. |
---|
| 226 | int startIndex = 0; |
---|
[1065] | 227 | |
---|
[1286] | 228 | // At most "fixThreshold" number of integer variables can be unsatisfied |
---|
| 229 | // for calling smallBranchAndBound(). |
---|
| 230 | // The PF Heuristic only fixes fixThreshold number of variables to |
---|
| 231 | // their integer values. Not more. Not less. The reason is to give |
---|
| 232 | // the smallBB some opportunity to find better solutions. If we fix |
---|
| 233 | // everything it might be too many (leading the heuristic to come up |
---|
| 234 | // with infeasibility rather than a useful result). |
---|
| 235 | // (This is an important paramater. And it is dynamically set.) |
---|
| 236 | double fixThreshold; |
---|
| 237 | /* |
---|
| 238 | if(numInt > 400) |
---|
| 239 | fixThreshold = 17*sqrt(numInt); |
---|
| 240 | if(numInt<=400 && numInt>100) |
---|
| 241 | fixThreshold = 5*sqrt(numInt); |
---|
| 242 | if(numInt<=100) |
---|
| 243 | fixThreshold = 4*sqrt(numInt); |
---|
| 244 | */ |
---|
| 245 | // Initialize fixThreshold based on the number of integer |
---|
| 246 | // variables |
---|
| 247 | if (numInt <= 100) |
---|
| 248 | fixThreshold = .35 * numInt; |
---|
| 249 | if (numInt > 100 && numInt < 1000) |
---|
| 250 | fixThreshold = .85 * numInt; |
---|
| 251 | if (numInt >= 1000) |
---|
| 252 | fixThreshold = .1 * numInt; |
---|
[1065] | 253 | |
---|
[1286] | 254 | // Whenever the dynamic system for changing fixThreshold |
---|
| 255 | // kicks in, it changes the parameter by the |
---|
| 256 | // fixThresholdChange amount. |
---|
| 257 | // (The 25% should be member data and tuned. Another paper!) |
---|
| 258 | double fixThresholdChange = 0.25 * fixThreshold; |
---|
[1065] | 259 | |
---|
[1286] | 260 | // maxNode is the maximum number of nodes we allow smallBB to |
---|
| 261 | // search. It's initialized to 400 and changed dynamically. |
---|
| 262 | // The 400 should be member data, if we become virtuous. |
---|
| 263 | int maxNode = 400; |
---|
[1065] | 264 | |
---|
[1286] | 265 | // We control the decision to change maxNode through the boolean |
---|
| 266 | // variable changeMaxNode. The boolean variable is initialized to |
---|
| 267 | // true and gets set to false under a condition (and is never true |
---|
| 268 | // again.) |
---|
| 269 | // It's flipped off and stays off (in the current incarnation of PF) |
---|
| 270 | bool changeMaxNode = 1; |
---|
[1065] | 271 | |
---|
[1286] | 272 | // The sumReturnCode is used for the dynamic system that sets |
---|
| 273 | // fixThreshold and changeMaxNode. |
---|
| 274 | // |
---|
| 275 | // We track what's happening in sumReturnCode. There are 8 switches. |
---|
| 276 | // The first 5 switches corresponds to a return code for smallBB. |
---|
| 277 | // |
---|
| 278 | // We want to know how many times we consecutively get the same |
---|
| 279 | // return code. |
---|
| 280 | // |
---|
| 281 | // If "good" return codes are happening often enough, we're happy. |
---|
| 282 | // |
---|
| 283 | // If a "bad" returncodes happen consecutively, we want to |
---|
| 284 | // change something. |
---|
| 285 | // |
---|
| 286 | // The switch 5 is the number of times PF didn't call smallBB |
---|
| 287 | // becuase the number of integer variables that took integer values |
---|
| 288 | // was less than fixThreshold. |
---|
| 289 | // |
---|
| 290 | // The swicth 6 was added for a brilliant idea...to be announced |
---|
| 291 | // later (another paper!) |
---|
| 292 | // |
---|
| 293 | // The switch 7 is the one that changes the max node. Read the |
---|
| 294 | // code. (Todo: Verbalize the brilliant idea for the masses.) |
---|
| 295 | // |
---|
| 296 | int sumReturnCode[8]; |
---|
| 297 | /* |
---|
| 298 | sumReturnCode[0] ~ -1 --> problem too big for smallBB |
---|
| 299 | sumReturnCode[1] ~ 0 --> smallBB not finshed and no soln |
---|
| 300 | sumReturnCode[2] ~ 1 --> smallBB not finshed and there is a soln |
---|
| 301 | sumReturnCode[3] ~ 2 --> smallBB finished and no soln |
---|
| 302 | sumReturnCode[4] ~ 3 --> smallBB finished and there is a soln |
---|
| 303 | sumReturnCode[5] ~ didn't call smallBranchAndBound too few to fix |
---|
| 304 | sumReturnCode[6] ~ didn't call smallBranchAndBound too many unsatisfied |
---|
| 305 | sumReturnCode[7] ~ the same as sumReturnCode[1] but becomes zero just if the returnCode is not 0 |
---|
| 306 | */ |
---|
[1065] | 307 | |
---|
[1286] | 308 | for (int i = 0; i < 8; i++) |
---|
| 309 | sumReturnCode[i] = 0; |
---|
| 310 | int * colIndex = new int[numCols]; |
---|
| 311 | for (int i = 0; i < numCols; i++) |
---|
| 312 | colIndex[i] = i; |
---|
| 313 | double cutoff = COIN_DBL_MAX; |
---|
| 314 | bool didMiniBB; |
---|
[1065] | 315 | |
---|
[1286] | 316 | // Main loop |
---|
| 317 | for (int i = 0; i < numRows; i++) { |
---|
| 318 | // track the number of mini-bb for the dynamic threshold setting |
---|
| 319 | didMiniBB = 0; |
---|
[1065] | 320 | |
---|
[1286] | 321 | for (int k = startIndex; k < solutions.numberSolutions; k++) |
---|
| 322 | //if the point has 0 unsatisfied variables; make sure it is |
---|
| 323 | //feasible. Check integer feasiblity and constraints. |
---|
| 324 | if (solutions.numberUnsatisfied[k] == 0) { |
---|
| 325 | double feasibility = 1; |
---|
| 326 | //check integer feasibility |
---|
| 327 | for (int icol = 0; icol < numCols; icol++) { |
---|
| 328 | double closest = floor(solutions.solution[k][icol] + 0.5); |
---|
| 329 | if (varClassInt[icol] && (fabs(solutions.solution[k][icol] - closest) > 1e-6)) { |
---|
| 330 | feasibility = 0; |
---|
| 331 | break; |
---|
| 332 | } |
---|
| 333 | } |
---|
| 334 | //check if the solution satisfies the constraints |
---|
| 335 | for (int irow = 0; irow < numRows; irow++) { |
---|
| 336 | double lhs = 0; |
---|
| 337 | for (int j = 0; j < numCols; j++) |
---|
| 338 | lhs += matrix[irow][j] * solutions.solution[k][j]; |
---|
| 339 | if (rowSense[irow] == 'L' && lhs > rhs[irow] + 1e-6) { |
---|
| 340 | feasibility = 0; |
---|
| 341 | break; |
---|
| 342 | } |
---|
| 343 | if (rowSense[irow] == 'G' && lhs < rhs[irow] - 1e-6) { |
---|
| 344 | feasibility = 0; |
---|
| 345 | break; |
---|
| 346 | } |
---|
| 347 | if (rowSense[irow] == 'E' && (lhs - rhs[irow] > 1e-6 || lhs - rhs[irow] < -1e-6)) { |
---|
| 348 | feasibility = 0; |
---|
| 349 | break; |
---|
| 350 | } |
---|
| 351 | } |
---|
[1065] | 352 | |
---|
[1286] | 353 | //if feasible, find the objective value and set the cutoff |
---|
| 354 | // for the smallBB and add a new constraint to the LP |
---|
| 355 | // (and update the best solution found so far for the |
---|
| 356 | // return arguments) |
---|
| 357 | if (feasibility) { |
---|
| 358 | double objectiveValue = 0; |
---|
| 359 | for (int j = 0; j < numCols; j++) |
---|
| 360 | objectiveValue += solutions.solution[k][j] * originalObjCoeff[j]; |
---|
| 361 | cutoff = objectiveValue; |
---|
| 362 | clpSolver->addRow(numCols, addRowIndex, originalObjCoeff, -COIN_DBL_MAX, cutoff); |
---|
[1065] | 363 | |
---|
[1286] | 364 | // Todo: pick up the best solution in the block (not |
---|
| 365 | // the last). |
---|
| 366 | solutionValue = objectiveValue; |
---|
| 367 | for (int m = 0; m < numCols; m++) |
---|
| 368 | betterSolution[m] = solutions.solution[k][m]; |
---|
| 369 | numFeasibles++; |
---|
| 370 | } |
---|
| 371 | } |
---|
[1065] | 372 | |
---|
[1286] | 373 | // Go through the block of solution and decide if to call smallBB |
---|
| 374 | for (int k = startIndex; k < solutions.numberSolutions; k++) { |
---|
| 375 | if (solutions.numberUnsatisfied[k] <= fixThreshold) { |
---|
| 376 | // get new copy |
---|
| 377 | OsiSolverInterface * newSolver; |
---|
| 378 | newSolver = new OsiClpSolverInterface(*clpSolver); |
---|
| 379 | newSolver->setObjSense(1); |
---|
| 380 | newSolver->setObjective(originalObjCoeff); |
---|
| 381 | int numberColumns = newSolver->getNumCols(); |
---|
| 382 | int numFixed = 0; |
---|
[1065] | 383 | |
---|
[1286] | 384 | // Fix the first fixThreshold number of integer vars |
---|
| 385 | // that are satisfied |
---|
| 386 | for (int iColumn = 0 ; iColumn < numberColumns ; iColumn++) { |
---|
| 387 | if (newSolver->isInteger(iColumn)) { |
---|
| 388 | double value = solutions.solution[k][iColumn]; |
---|
| 389 | double intValue = floor(value + 0.5); |
---|
| 390 | if (fabs(value - intValue) < 1.0e-5) { |
---|
| 391 | newSolver->setColLower(iColumn, intValue); |
---|
| 392 | newSolver->setColUpper(iColumn, intValue); |
---|
| 393 | numFixed++; |
---|
| 394 | if (numFixed > numInt - fixThreshold) |
---|
| 395 | break; |
---|
| 396 | } |
---|
| 397 | } |
---|
| 398 | } |
---|
[1641] | 399 | COIN_DETAIL_PRINT(printf("numFixed: %d\n", numFixed)); |
---|
| 400 | COIN_DETAIL_PRINT(printf("fixThreshold: %f\n", fixThreshold)); |
---|
| 401 | COIN_DETAIL_PRINT(printf("numInt: %d\n", numInt)); |
---|
[1286] | 402 | double *newSolution = new double[numCols]; |
---|
| 403 | double newSolutionValue; |
---|
[1065] | 404 | |
---|
[1286] | 405 | // Call smallBB on the modified problem |
---|
| 406 | int returnCode = smallBranchAndBound(newSolver, maxNode, newSolution, |
---|
| 407 | newSolutionValue, cutoff, "mini"); |
---|
| 408 | |
---|
| 409 | // If smallBB found a solution, update the better |
---|
| 410 | // solution and solutionValue (we gave smallBB our |
---|
| 411 | // cutoff, so it only finds improving solutions) |
---|
| 412 | if (returnCode == 1 || returnCode == 3) { |
---|
| 413 | numFeasibles ++; |
---|
| 414 | solutionValue = newSolutionValue; |
---|
| 415 | for (int m = 0; m < numCols; m++) |
---|
| 416 | betterSolution[m] = newSolution[m]; |
---|
[1641] | 417 | COIN_DETAIL_PRINT(printf("cutoff: %f\n", newSolutionValue)); |
---|
| 418 | COIN_DETAIL_PRINT(printf("time: %.2lf\n", CoinCpuTime() - start)); |
---|
[1286] | 419 | } |
---|
| 420 | didMiniBB = 1; |
---|
[1641] | 421 | COIN_DETAIL_PRINT(printf("returnCode: %d\n", returnCode)); |
---|
[1286] | 422 | |
---|
| 423 | //Update sumReturnCode array |
---|
| 424 | for (int iRC = 0; iRC < 6; iRC++) { |
---|
| 425 | if (iRC == returnCode + 1) |
---|
| 426 | sumReturnCode[iRC]++; |
---|
| 427 | else |
---|
| 428 | sumReturnCode[iRC] = 0; |
---|
| 429 | } |
---|
| 430 | if (returnCode != 0) |
---|
| 431 | sumReturnCode[7] = 0; |
---|
| 432 | else |
---|
| 433 | sumReturnCode[7]++; |
---|
| 434 | if (returnCode == 1 || returnCode == 3) { |
---|
| 435 | cutoff = newSolutionValue; |
---|
| 436 | clpSolver->addRow(numCols, addRowIndex, originalObjCoeff, -COIN_DBL_MAX, cutoff); |
---|
[1641] | 437 | COIN_DETAIL_PRINT(printf("******************\n\n*****************\n")); |
---|
[1286] | 438 | } |
---|
| 439 | break; |
---|
| 440 | } |
---|
| 441 | } |
---|
| 442 | |
---|
| 443 | if (!didMiniBB && solutions.numberSolutions - startIndex > 0) { |
---|
| 444 | sumReturnCode[5]++; |
---|
| 445 | for (int iRC = 0; iRC < 5; iRC++) |
---|
| 446 | sumReturnCode[iRC] = 0; |
---|
| 447 | } |
---|
| 448 | |
---|
| 449 | //Change "fixThreshold" if needed |
---|
| 450 | // using the data we've recorded in sumReturnCode |
---|
| 451 | if (sumReturnCode[1] >= 3) |
---|
| 452 | fixThreshold -= fixThresholdChange; |
---|
| 453 | if (sumReturnCode[7] >= 3 && changeMaxNode) { |
---|
| 454 | maxNode *= 5; |
---|
| 455 | changeMaxNode = 0; |
---|
| 456 | } |
---|
| 457 | if (sumReturnCode[3] >= 3 && fixThreshold < 0.95 * numInt) |
---|
| 458 | fixThreshold += fixThresholdChange; |
---|
| 459 | if (sumReturnCode[5] >= 4) |
---|
| 460 | fixThreshold += fixThresholdChange; |
---|
| 461 | if (sumReturnCode[0] > 3) |
---|
| 462 | fixThreshold -= fixThresholdChange; |
---|
| 463 | |
---|
| 464 | startIndex = solutions.numberSolutions; |
---|
| 465 | |
---|
| 466 | //Check if the maximum iterations limit is reached |
---|
| 467 | // rlh: Ask John how this is working with the change to trustedUserPtr. |
---|
| 468 | if (solutions.numberSolutions > 20000) |
---|
| 469 | break; |
---|
| 470 | |
---|
| 471 | // The first time in this loop PF solves orig LP. |
---|
| 472 | |
---|
| 473 | //Generate the random objective function |
---|
| 474 | randNum = rand() % 10 + 1; |
---|
| 475 | randNum = fmod(randNum, 2); |
---|
| 476 | for (int j = 0; j < numCols; j++) { |
---|
| 477 | if (randNum == 1) |
---|
| 478 | if (fabs(matrix[index[i]][j]) < 1e-6) |
---|
| 479 | newObj[j] = 0.1; |
---|
| 480 | else |
---|
| 481 | newObj[j] = matrix[index[i]][j] * 1.1; |
---|
| 482 | else if (fabs(matrix[index[i]][j]) < 1e-6) |
---|
| 483 | newObj[j] = -0.1; |
---|
| 484 | else |
---|
| 485 | newObj[j] = matrix[index[i]][j] * 0.9; |
---|
| 486 | } |
---|
| 487 | clpSolver->setObjective(newObj); |
---|
| 488 | if (rowSense[i] == 'L') |
---|
| 489 | clpSolver->setObjSense(-1); |
---|
| 490 | else |
---|
| 491 | // Todo #1: We don't need to solve the LPs to optimality. |
---|
| 492 | // We just need corner points. |
---|
| 493 | // There's a problem in stopping Clp that needs to be looked |
---|
| 494 | // into. So for now, we solve optimality. |
---|
| 495 | clpSolver->setObjSense(1); |
---|
| 496 | // simplex->setMaximumIterations(100); |
---|
| 497 | clpSolver->getModelPtr()->primal(1); |
---|
| 498 | // simplex->setMaximumIterations(100000); |
---|
[1641] | 499 | #ifdef COIN_DETAIL |
---|
[1286] | 500 | printf("cutoff: %f\n", cutoff); |
---|
| 501 | printf("time: %.2f\n", CoinCpuTime() - start); |
---|
| 502 | for (int iRC = 0; iRC < 8; iRC++) |
---|
| 503 | printf("%d ", sumReturnCode[iRC]); |
---|
| 504 | printf("\nfixThreshold: %f\n", fixThreshold); |
---|
| 505 | printf("numInt: %d\n", numInt); |
---|
| 506 | printf("\n---------------------------------------------------------------- %d\n", i); |
---|
[1641] | 507 | #endif |
---|
[1286] | 508 | |
---|
| 509 | //temp: |
---|
| 510 | if (i > 3) break; |
---|
| 511 | |
---|
[1065] | 512 | } |
---|
| 513 | |
---|
[1641] | 514 | COIN_DETAIL_PRINT(printf("Best Feasible Found: %f\n", cutoff)); |
---|
| 515 | COIN_DETAIL_PRINT(printf("Total time: %.2f\n", CoinCpuTime() - start)); |
---|
[1065] | 516 | |
---|
[1286] | 517 | if (numFeasibles == 0) { |
---|
| 518 | return 0; |
---|
| 519 | } |
---|
[1065] | 520 | |
---|
| 521 | |
---|
[1286] | 522 | |
---|
| 523 | // We found something better |
---|
| 524 | std::cout << "See you soon! You're leaving the Pivot-and-Fix Heuristic" << std::endl; |
---|
| 525 | std::cout << std::endl; |
---|
| 526 | |
---|
| 527 | return 1; |
---|
[1065] | 528 | #endif |
---|
| 529 | |
---|
[1286] | 530 | return 0; |
---|
[1065] | 531 | |
---|
[1055] | 532 | } |
---|
[1065] | 533 | |
---|
| 534 | |
---|
| 535 | |
---|
| 536 | |
---|
[1055] | 537 | // update model |
---|
[1271] | 538 | void CbcHeuristicPivotAndFix::setModel(CbcModel * ) |
---|
[1055] | 539 | { |
---|
[1286] | 540 | // probably same as resetModel |
---|
[1055] | 541 | } |
---|
| 542 | |
---|
[1286] | 543 | |
---|