Changeset 439 for branches/devel/Cbc/src/CbcHeuristicFPump.cpp
- Timestamp:
- Oct 8, 2006 7:33:47 PM (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/devel/Cbc/src/CbcHeuristicFPump.cpp
r435 r439 23 23 startTime_(0.0), 24 24 maximumTime_(0.0), 25 downValue_(0.5), 26 fakeCutoff_(COIN_DBL_MAX), 27 absoluteIncrement_(0.0), 28 relativeIncrement_(0.0), 25 29 maximumPasses_(100), 26 downValue_(0.5),30 maximumRetries_(1), 27 31 roundExpensive_(false) 28 32 { … … 36 40 startTime_(0.0), 37 41 maximumTime_(0.0), 42 downValue_(downValue), 43 fakeCutoff_(COIN_DBL_MAX), 44 absoluteIncrement_(0.0), 45 relativeIncrement_(0.0), 38 46 maximumPasses_(100), 39 downValue_(downValue),47 maximumRetries_(1), 40 48 roundExpensive_(roundExpensive) 41 49 { … … 65 73 else 66 74 fprintf(fp,"4 heuristicFPump.setMaximumPasses(%d);\n",maximumPasses_); 75 if (maximumRetries_!=other.maximumRetries_) 76 fprintf(fp,"3 heuristicFPump.setMaximumRetries(%d);\n",maximumRetries_); 77 else 78 fprintf(fp,"4 heuristicFPump.setMaximumRetries(%d);\n",maximumRetries_); 67 79 if (maximumTime_!=other.maximumTime_) 68 80 fprintf(fp,"3 heuristicFPump.setMaximumTime(%g);\n",maximumTime_); 69 81 else 70 82 fprintf(fp,"4 heuristicFPump.setMaximumTime(%g);\n",maximumTime_); 83 if (fakeCutoff_!=other.fakeCutoff_) 84 fprintf(fp,"3 heuristicFPump.setFakeCutoff(%g);\n",fakeCutoff_); 85 else 86 fprintf(fp,"4 heuristicFPump.setFakeCutoff(%g);\n",fakeCutoff_); 87 if (absoluteIncrement_!=other.absoluteIncrement_) 88 fprintf(fp,"3 heuristicFPump.setAbsoluteIncrement(%g);\n",absoluteIncrement_); 89 else 90 fprintf(fp,"4 heuristicFPump.setAbsoluteIncrement(%g);\n",absoluteIncrement_); 91 if (relativeIncrement_!=other.relativeIncrement_) 92 fprintf(fp,"3 heuristicFPump.setRelativeIncrement(%g);\n",relativeIncrement_); 93 else 94 fprintf(fp,"4 heuristicFPump.setRelativeIncrement(%g);\n",relativeIncrement_); 71 95 fprintf(fp,"3 cbcModel->addHeuristic(&heuristicFPump);\n"); 72 96 } … … 78 102 startTime_(rhs.startTime_), 79 103 maximumTime_(rhs.maximumTime_), 104 downValue_(rhs.downValue_), 105 fakeCutoff_(rhs.fakeCutoff_), 106 absoluteIncrement_(rhs.absoluteIncrement_), 107 relativeIncrement_(rhs.relativeIncrement_), 80 108 maximumPasses_(rhs.maximumPasses_), 81 downValue_(rhs.downValue_),109 maximumRetries_(rhs.maximumRetries_), 82 110 roundExpensive_(rhs.roundExpensive_) 83 111 { … … 109 137 // probably a good idea 110 138 if (model_->getSolutionCount()) return 0; 111 // Clone solver - otherwise annoys root node computations 112 OsiSolverInterface * solver = model_->solver()->clone(); 113 solver->setDblParam(OsiDualObjectiveLimit,1.0e50); 114 solver->resolve(); 115 const double * lower = solver->getColLower(); 116 const double * upper = solver->getColUpper(); 117 const double * solution = solver->getColSolution(); 118 double primalTolerance; 119 solver->getDblParam(OsiPrimalTolerance,primalTolerance); 120 139 // loop round doing repeated pumps 140 double cutoff; 141 model_->solver()->getDblParam(OsiDualObjectiveLimit,cutoff); 142 double direction = model_->solver()->getObjSense(); 143 cutoff *= direction; 144 cutoff = CoinMin(cutoff,solutionValue); 145 // space for rounded solution 121 146 int numberColumns = model_->getNumCols(); 147 double * newSolution = new double [numberColumns]; 148 double newSolutionValue=COIN_DBL_MAX; 149 bool solutionFound=false; 122 150 char * usedColumn = NULL; 123 151 double * lastSolution=NULL; … … 130 158 usedColumn = new char [numberColumns]; 131 159 memset(usedColumn,0,numberColumns); 132 lastSolution = CoinCopyOfArray(solver->getColSolution(),numberColumns); 160 model_->solver()->resolve(); 161 assert (model_->solver()->isProvenOptimal()); 162 lastSolution = CoinCopyOfArray(model_->solver()->getColSolution(),numberColumns); 133 163 } 134 int numberIntegers = model_->numberIntegers(); 135 const int * integerVariable = model_->integerVariable(); 136 137 // 1. initially check 0-1 138 int i,j; 139 int general=0; 140 for (i=0;i<numberIntegers;i++) { 141 int iColumn = integerVariable[i]; 164 int finalReturnCode=0; 165 int totalNumberPasses=0; 166 int numberTries=0; 167 while (true) { 168 int numberPasses=0; 169 numberTries++; 170 // Clone solver - otherwise annoys root node computations 171 OsiSolverInterface * solver = model_->solver()->clone(); 172 // if cutoff exists then add constraint 173 if (fabs(cutoff)<1.0e20&&(fakeCutoff_!=COIN_DBL_MAX||numberTries>1)) { 174 cutoff = CoinMin(cutoff,fakeCutoff_); 175 const double * objective = solver->getObjCoefficients(); 176 int numberColumns = solver->getNumCols(); 177 int * which = new int[numberColumns]; 178 double * els = new double[numberColumns]; 179 int nel=0; 180 for (int i=0;i<numberColumns;i++) { 181 double value = objective[i]; 182 if (value) { 183 which[nel]=i; 184 els[nel++] = direction*value; 185 } 186 } 187 solver->addRow(nel,which,els,-COIN_DBL_MAX,cutoff); 188 delete [] which; 189 delete [] els; 190 } 191 solver->setDblParam(OsiDualObjectiveLimit,1.0e50); 192 solver->resolve(); 193 // Solver may not be feasible 194 if (!solver->isProvenOptimal()) { 195 break; 196 } 197 const double * lower = solver->getColLower(); 198 const double * upper = solver->getColUpper(); 199 const double * solution = solver->getColSolution(); 200 if (lastSolution) 201 memcpy(lastSolution,solution,numberColumns*sizeof(double)); 202 double primalTolerance; 203 solver->getDblParam(OsiPrimalTolerance,primalTolerance); 204 205 int numberIntegers = model_->numberIntegers(); 206 const int * integerVariable = model_->integerVariable(); 207 208 // 1. initially check 0-1 209 int i,j; 210 int general=0; 211 for (i=0;i<numberIntegers;i++) { 212 int iColumn = integerVariable[i]; 142 213 #ifndef NDEBUG 143 const CbcObject * object = model_->object(i); 144 const CbcSimpleInteger * integerObject = 145 dynamic_cast<const CbcSimpleInteger *> (object); 146 assert(integerObject); 214 const OsiObject * object = model_->object(i); 215 const CbcSimpleInteger * integerObject = 216 dynamic_cast<const CbcSimpleInteger *> (object); 217 const OsiSimpleInteger * integerObject2 = 218 dynamic_cast<const OsiSimpleInteger *> (object); 219 assert(integerObject||integerObject2); 147 220 #endif 148 if (upper[iColumn]-lower[iColumn]>1.000001) { 149 general++; 150 break; 151 } 152 } 153 if (general*3>numberIntegers) { 154 delete solver; 155 return 0; 156 } 157 158 // 2. space for rounded solution 159 double * newSolution = new double [numberColumns]; 160 // space for last rounded solutions 221 if (upper[iColumn]-lower[iColumn]>1.000001) { 222 general++; 223 break; 224 } 225 } 226 if (general*3>numberIntegers) { 227 delete solver; 228 return 0; 229 } 230 231 // 2 space for last rounded solutions 161 232 #define NUMBER_OLD 4 162 double ** oldSolution = new double * [NUMBER_OLD]; 163 for (j=0;j<NUMBER_OLD;j++) { 164 oldSolution[j]= new double[numberColumns]; 165 for (i=0;i<numberColumns;i++) oldSolution[j][i]=-COIN_DBL_MAX; 166 } 167 168 // 3. Replace objective with an initial 0-valued objective 169 double * saveObjective = new double [numberColumns]; 170 memcpy(saveObjective,solver->getObjCoefficients(),numberColumns*sizeof(double)); 171 for (i=0;i<numberColumns;i++) { 172 solver->setObjCoeff(i,0.0); 173 } 174 bool finished=false; 175 double direction = solver->getObjSense(); 176 int returnCode=0; 177 bool takeHint; 178 OsiHintStrength strength; 179 solver->getHintParam(OsiDoDualInResolve,takeHint,strength); 180 solver->setHintParam(OsiDoDualInResolve,false); 181 solver->messageHandler()->setLogLevel(0); 182 183 // 4. Save objective offset so we can see progress 184 double saveOffset; 185 solver->getDblParam(OsiObjOffset,saveOffset); 186 187 // 5. MAIN WHILE LOOP 188 int numberPasses=0; 189 double newSolutionValue=COIN_DBL_MAX; 190 bool newLineNeeded=false; 191 while (!finished) { 192 returnCode=0; 193 // see what changed 194 if (usedColumn) { 195 const double * solution = solver->getColSolution(); 196 for (i=0;i<numberColumns;i++) { 197 if (fabs(solution[i]-lastSolution[i])>1.0e-8) 198 usedColumn[i]=1; 199 lastSolution[i]=solution[i]; 200 } 201 } 202 if (numberPasses>=maximumPasses_) { 203 break; 204 } 205 if (maximumTime_>0.0&&CoinCpuTime()>=startTime_+maximumTime_) break; 206 numberPasses++; 207 memcpy(newSolution,solution,numberColumns*sizeof(double)); 208 int flip; 209 returnCode = rounds(newSolution,saveObjective,roundExpensive_,downValue_,&flip); 210 if (returnCode) { 211 // SOLUTION IS INTEGER 212 // Put back correct objective 213 for (i=0;i<numberColumns;i++) 214 solver->setObjCoeff(i,saveObjective[i]); 215 // solution - but may not be better 216 // Compute using dot product 217 solver->setDblParam(OsiObjOffset,saveOffset); 218 newSolutionValue = direction*solver->OsiSolverInterface::getObjValue(); 219 if (model_->logLevel()) 220 printf(" - solution found\n"); 221 newLineNeeded=false; 222 if (newSolutionValue<solutionValue) { 223 if (general) { 224 int numberLeft=0; 225 for (i=0;i<numberIntegers;i++) { 226 int iColumn = integerVariable[i]; 227 double value = floor(newSolution[iColumn]+0.5); 228 if(solver->isBinary(iColumn)) { 229 solver->setColLower(iColumn,value); 230 solver->setColUpper(iColumn,value); 231 } else { 232 if (fabs(value-newSolution[iColumn])>1.0e-7) 233 numberLeft++; 234 } 235 } 236 if (numberLeft) { 237 returnCode = smallBranchAndBound(solver,200,newSolution,newSolutionValue, 238 solutionValue,"CbcHeuristicFpump"); 239 } 240 } 241 if (returnCode) { 242 memcpy(betterSolution,newSolution,numberColumns*sizeof(double)); 243 solutionValue=newSolutionValue; 244 } 233 double ** oldSolution = new double * [NUMBER_OLD]; 234 for (j=0;j<NUMBER_OLD;j++) { 235 oldSolution[j]= new double[numberColumns]; 236 for (i=0;i<numberColumns;i++) oldSolution[j][i]=-COIN_DBL_MAX; 237 } 238 239 // 3. Replace objective with an initial 0-valued objective 240 double * saveObjective = new double [numberColumns]; 241 memcpy(saveObjective,solver->getObjCoefficients(),numberColumns*sizeof(double)); 242 for (i=0;i<numberColumns;i++) { 243 solver->setObjCoeff(i,0.0); 244 } 245 bool finished=false; 246 double direction = solver->getObjSense(); 247 int returnCode=0; 248 bool takeHint; 249 OsiHintStrength strength; 250 solver->getHintParam(OsiDoDualInResolve,takeHint,strength); 251 solver->setHintParam(OsiDoDualInResolve,false); 252 solver->messageHandler()->setLogLevel(0); 253 254 // 4. Save objective offset so we can see progress 255 double saveOffset; 256 solver->getDblParam(OsiObjOffset,saveOffset); 257 258 // 5. MAIN WHILE LOOP 259 bool newLineNeeded=false; 260 while (!finished) { 261 returnCode=0; 262 // see what changed 263 if (usedColumn) { 264 for (i=0;i<numberColumns;i++) { 265 if (fabs(solution[i]-lastSolution[i])>1.0e-8) 266 usedColumn[i]=1; 267 lastSolution[i]=solution[i]; 268 } 269 } 270 if (numberPasses>=maximumPasses_) { 271 break; 272 } 273 if (maximumTime_>0.0&&CoinCpuTime()>=startTime_+maximumTime_) break; 274 numberPasses++; 275 memcpy(newSolution,solution,numberColumns*sizeof(double)); 276 int flip; 277 returnCode = rounds(newSolution,saveObjective,roundExpensive_,downValue_,&flip); 278 if (returnCode) { 279 // SOLUTION IS INTEGER 280 // Put back correct objective 281 for (i=0;i<numberColumns;i++) 282 solver->setObjCoeff(i,saveObjective[i]); 283 284 // solution - but may not be better 285 // Compute using dot product 286 solver->setDblParam(OsiObjOffset,saveOffset); 287 newSolutionValue = -saveOffset; 288 for ( i=0 ; i<numberColumns ; i++ ) 289 newSolutionValue += saveObjective[i]*newSolution[i]; 290 newSolutionValue *= direction; 291 if (model_->logLevel()) 292 printf(" - solution found of %g",newSolutionValue); 293 newLineNeeded=false; 294 if (newSolutionValue<solutionValue) { 295 double saveValue = newSolutionValue; 296 if (general) { 297 int numberLeft=0; 298 for (i=0;i<numberIntegers;i++) { 299 int iColumn = integerVariable[i]; 300 double value = floor(newSolution[iColumn]+0.5); 301 if(solver->isBinary(iColumn)) { 302 solver->setColLower(iColumn,value); 303 solver->setColUpper(iColumn,value); 304 } else { 305 if (fabs(value-newSolution[iColumn])>1.0e-7) 306 numberLeft++; 307 } 308 } 309 if (numberLeft) { 310 returnCode = smallBranchAndBound(solver,200,newSolution,newSolutionValue, 311 solutionValue,"CbcHeuristicFpump"); 312 } 313 } 314 if (returnCode) { 315 memcpy(betterSolution,newSolution,numberColumns*sizeof(double)); 316 solutionValue=newSolutionValue; 317 solutionFound=true; 318 if (general&&saveValue!=newSolutionValue) 319 printf(" - cleaned solution of %g\n",solutionValue); 320 else 321 printf("\n"); 322 } else { 323 if (model_->logLevel()) 324 printf(" - not good enough after small branch and bound\n"); 325 } 326 } else { 327 if (model_->logLevel()) 328 printf(" - not good enough\n"); 329 newLineNeeded=false; 330 returnCode=0; 331 } 332 break; 245 333 } else { 246 returnCode=0; 247 } 248 break; 249 } else { 250 // SOLUTION IS not INTEGER 251 // 1. check for loop 252 bool matched; 253 for (int k = NUMBER_OLD-1; k > 0; k--) { 334 // SOLUTION IS not INTEGER 335 // 1. check for loop 336 bool matched; 337 for (int k = NUMBER_OLD-1; k > 0; k--) { 254 338 double * b = oldSolution[k]; 255 339 matched = true; 256 340 for (i = 0; i <numberIntegers; i++) { 257 258 259 260 261 262 263 341 int iColumn = integerVariable[i]; 342 if(!solver->isBinary(iColumn)) 343 continue; 344 if (newSolution[iColumn]!=b[iColumn]) { 345 matched=false; 346 break; 347 } 264 348 } 265 349 if (matched) break; 266 } 267 if (matched || numberPasses%100 == 0) { 268 // perturbation 269 if (model_->logLevel()) 270 printf("Perturbation applied"); 271 newLineNeeded=true; 272 for (i=0;i<numberIntegers;i++) { 273 int iColumn = integerVariable[i]; 274 if(!solver->isBinary(iColumn)) 275 continue; 276 double value = max(0.0,CoinDrand48()-0.3); 277 double difference = fabs(solution[iColumn]-newSolution[iColumn]); 278 if (difference+value>0.5) { 279 if (newSolution[iColumn]<lower[iColumn]+primalTolerance) newSolution[iColumn] += 1.0; 280 else if (newSolution[iColumn]>upper[iColumn]-primalTolerance) newSolution[iColumn] -= 1.0; 281 else abort(); 282 } 283 } 284 } else { 285 for (j=NUMBER_OLD-1;j>0;j--) { 286 for (i = 0; i < numberColumns; i++) oldSolution[j][i]=oldSolution[j-1][i]; 287 } 288 for (j = 0; j < numberColumns; j++) oldSolution[0][j] = newSolution[j]; 289 } 290 291 // 2. update the objective function based on the new rounded solution 292 double offset=0.0; 350 } 351 if (matched || numberPasses%100 == 0) { 352 // perturbation 353 if (model_->logLevel()) 354 printf("Perturbation applied"); 355 newLineNeeded=true; 356 for (i=0;i<numberIntegers;i++) { 357 int iColumn = integerVariable[i]; 358 if(!solver->isBinary(iColumn)) 359 continue; 360 double value = max(0.0,CoinDrand48()-0.3); 361 double difference = fabs(solution[iColumn]-newSolution[iColumn]); 362 if (difference+value>0.5) { 363 if (newSolution[iColumn]<lower[iColumn]+primalTolerance) newSolution[iColumn] += 1.0; 364 else if (newSolution[iColumn]>upper[iColumn]-primalTolerance) newSolution[iColumn] -= 1.0; 365 else abort(); 366 } 367 } 368 } else { 369 for (j=NUMBER_OLD-1;j>0;j--) { 370 for (i = 0; i < numberColumns; i++) oldSolution[j][i]=oldSolution[j-1][i]; 371 } 372 for (j = 0; j < numberColumns; j++) oldSolution[0][j] = newSolution[j]; 373 } 374 375 // 2. update the objective function based on the new rounded solution 376 double offset=0.0; 377 for (i=0;i<numberIntegers;i++) { 378 int iColumn = integerVariable[i]; 379 if(!solver->isBinary(iColumn)) 380 continue; 381 double costValue = 1.0; 382 // deal with fixed variables (i.e., upper=lower) 383 if (fabs(lower[iColumn]-upper[iColumn]) < primalTolerance) { 384 if (lower[iColumn] > 1. - primalTolerance) solver->setObjCoeff(iColumn,-costValue); 385 else solver->setObjCoeff(iColumn,costValue); 386 continue; 387 } 388 if (newSolution[iColumn]<lower[iColumn]+primalTolerance) { 389 solver->setObjCoeff(iColumn,costValue); 390 } else { 391 if (newSolution[iColumn]>upper[iColumn]-primalTolerance) { 392 solver->setObjCoeff(iColumn,-costValue); 393 } else { 394 abort(); 395 } 396 } 397 offset += costValue*newSolution[iColumn]; 398 } 399 solver->setDblParam(OsiObjOffset,-offset); 400 if (!general&false) { 401 // Solve in two goes - first keep satisfied ones fixed 402 double * saveLower = new double [numberIntegers]; 403 double * saveUpper = new double [numberIntegers]; 404 for (i=0;i<numberIntegers;i++) { 405 int iColumn = integerVariable[i]; 406 saveLower[i]=COIN_DBL_MAX; 407 saveUpper[i]=-COIN_DBL_MAX; 408 if (solution[iColumn]<lower[iColumn]+primalTolerance) { 409 saveUpper[i]=upper[iColumn]; 410 solver->setColUpper(iColumn,lower[iColumn]); 411 } else if (solution[iColumn]>upper[iColumn]-primalTolerance) { 412 saveLower[i]=lower[iColumn]; 413 solver->setColLower(iColumn,upper[iColumn]); 414 } 415 } 416 solver->resolve(); 417 assert (solver->isProvenOptimal()); 418 for (i=0;i<numberIntegers;i++) { 419 int iColumn = integerVariable[i]; 420 if (saveLower[i]!=COIN_DBL_MAX) 421 solver->setColLower(iColumn,saveLower[i]); 422 if (saveUpper[i]!=-COIN_DBL_MAX) 423 solver->setColUpper(iColumn,saveUpper[i]); 424 saveUpper[i]=-COIN_DBL_MAX; 425 } 426 memcpy(newSolution,solution,numberColumns*sizeof(double)); 427 int flip; 428 returnCode = rounds(newSolution,saveObjective,roundExpensive_,downValue_,&flip); 429 if (returnCode) { 430 // solution - but may not be better 431 // Compute using dot product 432 double newSolutionValue = -saveOffset; 433 for ( i=0 ; i<numberColumns ; i++ ) 434 newSolutionValue += saveObjective[i]*newSolution[i]; 435 newSolutionValue *= direction; 436 if (model_->logLevel()) 437 printf(" - intermediate solution found of %g",newSolutionValue); 438 if (newSolutionValue<solutionValue) { 439 memcpy(betterSolution,newSolution,numberColumns*sizeof(double)); 440 solutionValue=newSolutionValue; 441 solutionFound=true; 442 } else { 443 returnCode=0; 444 } 445 } 446 } 447 solver->resolve(); 448 assert (solver->isProvenOptimal()); 449 if (model_->logLevel()) 450 printf("\npass %3d: obj. %10.5f --> ", numberPasses+totalNumberPasses,solver->getObjValue()); 451 newLineNeeded=true; 452 453 } 454 } // END WHILE 455 if (newLineNeeded&&model_->logLevel()) 456 printf(" - no solution found\n"); 457 delete solver; 458 for ( j=0;j<NUMBER_OLD;j++) 459 delete [] oldSolution[j]; 460 delete [] oldSolution; 461 delete [] saveObjective; 462 if (usedColumn) { 463 OsiSolverInterface * newSolver = model_->continuousSolver()->clone(); 464 const double * colLower = newSolver->getColLower(); 465 const double * colUpper = newSolver->getColUpper(); 466 int i; 467 int nFix=0; 468 int nFixI=0; 469 int nFixC=0; 293 470 for (i=0;i<numberIntegers;i++) { 294 int iColumn = integerVariable[i]; 295 if(!solver->isBinary(iColumn)) 296 continue; 297 double costValue = 1.0; 298 // deal with fixed variables (i.e., upper=lower) 299 if (fabs(lower[iColumn]-upper[iColumn]) < primalTolerance) { 300 if (lower[iColumn] > 1. - primalTolerance) solver->setObjCoeff(iColumn,-costValue); 301 else solver->setObjCoeff(iColumn,costValue); 302 continue; 303 } 304 if (newSolution[iColumn]<lower[iColumn]+primalTolerance) { 305 solver->setObjCoeff(iColumn,costValue); 306 } else { 307 if (newSolution[iColumn]>upper[iColumn]-primalTolerance) { 308 solver->setObjCoeff(iColumn,-costValue); 309 } else { 310 abort(); 311 } 312 } 313 offset += costValue*newSolution[iColumn]; 314 } 315 solver->setDblParam(OsiObjOffset,-offset); 316 solver->resolve(); 317 if (model_->logLevel()) 318 printf("\npass %3d: obj. %10.5f --> ", numberPasses,solver->getObjValue()); 319 newLineNeeded=true; 320 321 } 322 } // END WHILE 323 if (newLineNeeded&&model_->logLevel()) 324 printf(" - no solution found\n"); 325 delete solver; 326 for ( j=0;j<NUMBER_OLD;j++) 327 delete [] oldSolution[j]; 328 delete [] oldSolution; 329 delete [] saveObjective; 330 if (usedColumn) { 331 OsiSolverInterface * newSolver = model_->continuousSolver()->clone(); 332 const double * colLower = newSolver->getColLower(); 333 const double * colUpper = newSolver->getColUpper(); 334 int i; 335 int nFix=0; 336 int nFixI=0; 337 int nFixC=0; 338 for (i=0;i<numberIntegers;i++) { 339 int iColumn=integerVariable[i]; 340 const CbcObject * object = model_->object(i); 341 const CbcSimpleInteger * integerObject = 342 dynamic_cast<const CbcSimpleInteger *> (object); 343 assert(integerObject); 344 // get original bounds 345 double originalLower = integerObject->originalLowerBound(); 346 assert(colLower[iColumn]==originalLower); 347 newSolver->setColLower(iColumn,CoinMax(colLower[iColumn],originalLower)); 348 double originalUpper = integerObject->originalUpperBound(); 349 assert(colUpper[iColumn]==originalUpper); 350 newSolver->setColUpper(iColumn,CoinMin(colUpper[iColumn],originalUpper)); 351 if (!usedColumn[iColumn]) { 352 double value=lastSolution[iColumn]; 353 double nearest = floor(value+0.5); 354 if (fabs(value-nearest)<1.0e-7) { 355 if (nearest==colLower[iColumn]) { 356 newSolver->setColUpper(iColumn,colLower[iColumn]); 357 nFix++; 358 } else if (nearest==colUpper[iColumn]) { 359 newSolver->setColLower(iColumn,colUpper[iColumn]); 360 nFix++; 361 } else if (fixInternal) { 362 newSolver->setColLower(iColumn,nearest); 363 newSolver->setColUpper(iColumn,nearest); 364 nFix++; 365 nFixI++; 366 } 367 } 368 } 369 } 370 if (fixContinuous) { 371 for (int iColumn=0;iColumn<numberColumns;iColumn++) { 372 if (!newSolver->isInteger(iColumn)&&!usedColumn[iColumn]) { 471 int iColumn=integerVariable[i]; 472 const OsiObject * object = model_->object(i); 473 double originalLower; 474 double originalUpper; 475 getIntegerInformation( object,originalLower, originalUpper); 476 assert(colLower[iColumn]==originalLower); 477 newSolver->setColLower(iColumn,CoinMax(colLower[iColumn],originalLower)); 478 assert(colUpper[iColumn]==originalUpper); 479 newSolver->setColUpper(iColumn,CoinMin(colUpper[iColumn],originalUpper)); 480 if (!usedColumn[iColumn]) { 373 481 double value=lastSolution[iColumn]; 374 if (value<colLower[iColumn]+1.0e-8) { 375 newSolver->setColUpper(iColumn,colLower[iColumn]); 376 nFix++; 377 nFixC++; 378 } else if (value>colUpper[iColumn]-1.0e-8) { 379 newSolver->setColLower(iColumn,colUpper[iColumn]); 380 nFix++; 381 nFixC++; 382 } 383 } 384 } 385 } 386 newSolver->initialSolve(); 387 assert (newSolver->isProvenOptimal()); 388 printf("%d integers at bound fixed, %d internal and %d continuous\n", 389 nFix,nFixI,nFixC); 390 double saveValue = newSolutionValue; 391 returnCode = smallBranchAndBound(newSolver,200,newSolution,newSolutionValue, 392 newSolutionValue,"CbcHeuristicLocalAfterFPump"); 393 if (returnCode) { 394 printf("old sol of %g new of %g\n",saveValue,newSolutionValue); 395 memcpy(betterSolution,newSolution,numberColumns*sizeof(double)); 396 solutionValue=newSolutionValue; 397 } 398 delete newSolver; 399 delete [] usedColumn; 400 delete [] lastSolution; 482 double nearest = floor(value+0.5); 483 if (fabs(value-nearest)<1.0e-7) { 484 if (nearest==colLower[iColumn]) { 485 newSolver->setColUpper(iColumn,colLower[iColumn]); 486 nFix++; 487 } else if (nearest==colUpper[iColumn]) { 488 newSolver->setColLower(iColumn,colUpper[iColumn]); 489 nFix++; 490 } else if (fixInternal) { 491 newSolver->setColLower(iColumn,nearest); 492 newSolver->setColUpper(iColumn,nearest); 493 nFix++; 494 nFixI++; 495 } 496 } 497 } 498 } 499 if (fixContinuous) { 500 for (int iColumn=0;iColumn<numberColumns;iColumn++) { 501 if (!newSolver->isInteger(iColumn)&&!usedColumn[iColumn]) { 502 double value=lastSolution[iColumn]; 503 if (value<colLower[iColumn]+1.0e-8) { 504 newSolver->setColUpper(iColumn,colLower[iColumn]); 505 nFixC++; 506 } else if (value>colUpper[iColumn]-1.0e-8) { 507 newSolver->setColLower(iColumn,colUpper[iColumn]); 508 nFixC++; 509 } 510 } 511 } 512 } 513 newSolver->initialSolve(); 514 if (!newSolver->isProvenOptimal()) { 515 newSolver->writeMps("bad.mps"); 516 assert (newSolver->isProvenOptimal()); 517 break; 518 } 519 printf("%d integers at bound fixed (of which %d were internal) and %d continuous\n", 520 nFix,nFixI,nFixC); 521 double saveValue = newSolutionValue; 522 returnCode = smallBranchAndBound(newSolver,200,newSolution,newSolutionValue, 523 newSolutionValue,"CbcHeuristicLocalAfterFPump"); 524 if (returnCode) { 525 printf("old sol of %g new of %g\n",saveValue,newSolutionValue); 526 memcpy(betterSolution,newSolution,numberColumns*sizeof(double)); 527 solutionValue=newSolutionValue; 528 solutionFound=true; 529 } 530 delete newSolver; 531 } 532 if (solutionFound) finalReturnCode=1; 533 cutoff = CoinMin(cutoff,solutionValue); 534 if (numberTries>=CoinAbs(maximumRetries_)||!solutionFound) { 535 break; 536 } else if (absoluteIncrement_>0.0||relativeIncrement_>0.0) { 537 solutionFound=false; 538 double gap = relativeIncrement_*fabs(solutionValue); 539 cutoff -= CoinMax(CoinMax(gap,absoluteIncrement_),model_->getCutoffIncrement()); 540 printf("round again with cutoff of %g\n",cutoff); 541 if (maximumRetries_<0) 542 memset(usedColumn,0,numberColumns); 543 totalNumberPasses += numberPasses; 544 } else { 545 break; 546 } 401 547 } 548 delete [] usedColumn; 549 delete [] lastSolution; 402 550 delete [] newSolution; 403 return returnCode;551 return finalReturnCode; 404 552 } 405 553 … … 449 597 continue; 450 598 #ifndef NDEBUG 451 const CbcObject * object = model_->object(i); 452 const CbcSimpleInteger * integerObject = 453 dynamic_cast<const CbcSimpleInteger *> (object); 454 assert(integerObject); 599 const OsiObject * object = model_->object(i); 600 const CbcSimpleInteger * integerObject = 601 dynamic_cast<const CbcSimpleInteger *> (object); 602 const OsiSimpleInteger * integerObject2 = 603 dynamic_cast<const OsiSimpleInteger *> (object); 604 assert(integerObject||integerObject2); 455 605 #endif 456 606 double value=solution[iColumn];
Note: See TracChangeset
for help on using the changeset viewer.