389 | | int maxParam = 200; |
390 | | CbcOrClpParam * parameters = new CbcOrClpParam [maxParam]; |
391 | | numberParameters_=0 ; |
392 | | establishParams(numberParameters_,parameters) ; |
393 | | assert (numberParameters_<=maxParam); |
394 | | parameters_ = new CbcOrClpParam [numberParameters_]; |
395 | | int i; |
396 | | for (i=0;i<numberParameters_;i++) |
397 | | parameters_[i]=parameters[i]; |
398 | | delete [] parameters; |
399 | | const char dirsep = CoinFindDirSeparator(); |
400 | | std::string directory; |
401 | | std::string dirSample; |
402 | | std::string dirNetlib; |
403 | | std::string dirMiplib; |
404 | | if (dirsep == '/') { |
405 | | directory = "./"; |
406 | | dirSample = "../../Data/Sample/"; |
407 | | dirNetlib = "../../Data/Netlib/"; |
408 | | dirMiplib = "../../Data/miplib3/"; |
409 | | } else { |
410 | | directory = ".\\"; |
411 | | dirSample = "..\\..\\Data\\Sample\\"; |
412 | | dirNetlib = "..\\..\\Data\\Netlib\\"; |
413 | | dirMiplib = "..\\..\\Data\\miplib3\\"; |
414 | | } |
415 | | std::string defaultDirectory = directory; |
416 | | std::string importFile =""; |
417 | | std::string exportFile ="default.mps"; |
418 | | std::string importBasisFile =""; |
419 | | std::string importPriorityFile =""; |
420 | | std::string debugFile=""; |
421 | | std::string printMask=""; |
422 | | std::string exportBasisFile ="default.bas"; |
423 | | std::string saveFile ="default.prob"; |
424 | | std::string restoreFile ="default.prob"; |
425 | | std::string solutionFile ="stdout"; |
426 | | std::string solutionSaveFile ="solution.file"; |
427 | | int doIdiot=-1; |
428 | | int outputFormat=2; |
429 | | int substitution=3; |
430 | | int dualize=3; |
431 | | int preSolve=5; |
432 | | int doSprint=-1; |
433 | | int testOsiParameters=-1; |
434 | | int createSolver=0; |
435 | | ClpSimplex * lpSolver; |
436 | | OsiClpSolverInterface * clpSolver; |
437 | | if (model_.solver()) { |
438 | | clpSolver = dynamic_cast<OsiClpSolverInterface *> (model_.solver()); |
439 | | assert (clpSolver); |
440 | | lpSolver = clpSolver->getModelPtr(); |
441 | | assert (lpSolver); |
442 | | } else { |
443 | | lpSolver = new ClpSimplex(); |
444 | | clpSolver = new OsiClpSolverInterface(lpSolver,true); |
445 | | createSolver =1 ; |
446 | | } |
447 | | parameters_[whichParam(BASISIN,numberParameters_,parameters_)].setStringValue(importBasisFile); |
448 | | parameters_[whichParam(PRIORITYIN,numberParameters_,parameters_)].setStringValue(importPriorityFile); |
449 | | parameters_[whichParam(BASISOUT,numberParameters_,parameters_)].setStringValue(exportBasisFile); |
450 | | parameters_[whichParam(DEBUG,numberParameters_,parameters_)].setStringValue(debugFile); |
451 | | parameters_[whichParam(PRINTMASK,numberParameters_,parameters_)].setStringValue(printMask); |
452 | | parameters_[whichParam(DIRECTORY,numberParameters_,parameters_)].setStringValue(directory); |
453 | | parameters_[whichParam(DIRSAMPLE,numberParameters_,parameters_)].setStringValue(dirSample); |
454 | | parameters_[whichParam(DIRNETLIB,numberParameters_,parameters_)].setStringValue(dirNetlib); |
455 | | parameters_[whichParam(DIRMIPLIB,numberParameters_,parameters_)].setStringValue(dirMiplib); |
456 | | parameters_[whichParam(DUALBOUND,numberParameters_,parameters_)].setDoubleValue(lpSolver->dualBound()); |
457 | | parameters_[whichParam(DUALTOLERANCE,numberParameters_,parameters_)].setDoubleValue(lpSolver->dualTolerance()); |
458 | | parameters_[whichParam(EXPORT,numberParameters_,parameters_)].setStringValue(exportFile); |
459 | | parameters_[whichParam(IDIOT,numberParameters_,parameters_)].setIntValue(doIdiot); |
460 | | parameters_[whichParam(IMPORT,numberParameters_,parameters_)].setStringValue(importFile); |
461 | | parameters_[whichParam(PRESOLVETOLERANCE,numberParameters_,parameters_)].setDoubleValue(1.0e-8); |
462 | | int iParam = whichParam(SOLVERLOGLEVEL,numberParameters_,parameters_); |
463 | | int value=1; |
464 | | clpSolver->messageHandler()->setLogLevel(1) ; |
465 | | lpSolver->setLogLevel(1); |
466 | | parameters_[iParam].setIntValue(value); |
467 | | iParam = whichParam(LOGLEVEL,numberParameters_,parameters_); |
468 | | model_.messageHandler()->setLogLevel(value); |
469 | | parameters_[iParam].setIntValue(value); |
470 | | parameters_[whichParam(MAXFACTOR,numberParameters_,parameters_)].setIntValue(lpSolver->factorizationFrequency()); |
471 | | parameters_[whichParam(MAXITERATION,numberParameters_,parameters_)].setIntValue(lpSolver->maximumIterations()); |
472 | | parameters_[whichParam(OUTPUTFORMAT,numberParameters_,parameters_)].setIntValue(outputFormat); |
473 | | parameters_[whichParam(PRESOLVEPASS,numberParameters_,parameters_)].setIntValue(preSolve); |
474 | | parameters_[whichParam(PERTVALUE,numberParameters_,parameters_)].setIntValue(lpSolver->perturbation()); |
475 | | parameters_[whichParam(PRIMALTOLERANCE,numberParameters_,parameters_)].setDoubleValue(lpSolver->primalTolerance()); |
476 | | parameters_[whichParam(PRIMALWEIGHT,numberParameters_,parameters_)].setDoubleValue(lpSolver->infeasibilityCost()); |
477 | | parameters_[whichParam(RESTORE,numberParameters_,parameters_)].setStringValue(restoreFile); |
478 | | parameters_[whichParam(SAVE,numberParameters_,parameters_)].setStringValue(saveFile); |
479 | | //parameters_[whichParam(TIMELIMIT,numberParameters_,parameters_)].setDoubleValue(1.0e8); |
480 | | parameters_[whichParam(TIMELIMIT_BAB,numberParameters_,parameters_)].setDoubleValue(1.0e8); |
481 | | parameters_[whichParam(SOLUTION,numberParameters_,parameters_)].setStringValue(solutionFile); |
482 | | parameters_[whichParam(SAVESOL,numberParameters_,parameters_)].setStringValue(solutionSaveFile); |
483 | | parameters_[whichParam(SPRINT,numberParameters_,parameters_)].setIntValue(doSprint); |
484 | | parameters_[whichParam(SUBSTITUTION,numberParameters_,parameters_)].setIntValue(substitution); |
485 | | parameters_[whichParam(DUALIZE,numberParameters_,parameters_)].setIntValue(dualize); |
486 | | parameters_[whichParam(NUMBERBEFORE,numberParameters_,parameters_)].setIntValue(model_.numberBeforeTrust()); |
487 | | parameters_[whichParam(MAXNODES,numberParameters_,parameters_)].setIntValue(model_.getMaximumNodes()); |
488 | | parameters_[whichParam(STRONGBRANCHING,numberParameters_,parameters_)].setIntValue(model_.numberStrong()); |
489 | | parameters_[whichParam(INFEASIBILITYWEIGHT,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcInfeasibilityWeight)); |
490 | | parameters_[whichParam(INTEGERTOLERANCE,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcIntegerTolerance)); |
491 | | parameters_[whichParam(INCREMENT,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcCutoffIncrement)); |
492 | | parameters_[whichParam(TESTOSI,numberParameters_,parameters_)].setIntValue(testOsiParameters); |
493 | | parameters_[whichParam(FPUMPTUNE,numberParameters_,parameters_)].setIntValue(1003); |
494 | | initialPumpTune=1003; |
| 389 | int maxParam = 200; |
| 390 | CbcOrClpParam * parameters = new CbcOrClpParam [maxParam]; |
| 391 | numberParameters_ = 0 ; |
| 392 | establishParams(numberParameters_, parameters) ; |
| 393 | assert (numberParameters_ <= maxParam); |
| 394 | parameters_ = new CbcOrClpParam [numberParameters_]; |
| 395 | int i; |
| 396 | for (i = 0; i < numberParameters_; i++) |
| 397 | parameters_[i] = parameters[i]; |
| 398 | delete [] parameters; |
| 399 | const char dirsep = CoinFindDirSeparator(); |
| 400 | std::string directory; |
| 401 | std::string dirSample; |
| 402 | std::string dirNetlib; |
| 403 | std::string dirMiplib; |
| 404 | if (dirsep == '/') { |
| 405 | directory = "./"; |
| 406 | dirSample = "../../Data/Sample/"; |
| 407 | dirNetlib = "../../Data/Netlib/"; |
| 408 | dirMiplib = "../../Data/miplib3/"; |
| 409 | } else { |
| 410 | directory = ".\\"; |
| 411 | dirSample = "..\\..\\Data\\Sample\\"; |
| 412 | dirNetlib = "..\\..\\Data\\Netlib\\"; |
| 413 | dirMiplib = "..\\..\\Data\\miplib3\\"; |
| 414 | } |
| 415 | std::string defaultDirectory = directory; |
| 416 | std::string importFile = ""; |
| 417 | std::string exportFile = "default.mps"; |
| 418 | std::string importBasisFile = ""; |
| 419 | std::string importPriorityFile = ""; |
| 420 | std::string debugFile = ""; |
| 421 | std::string printMask = ""; |
| 422 | std::string exportBasisFile = "default.bas"; |
| 423 | std::string saveFile = "default.prob"; |
| 424 | std::string restoreFile = "default.prob"; |
| 425 | std::string solutionFile = "stdout"; |
| 426 | std::string solutionSaveFile = "solution.file"; |
| 427 | int doIdiot = -1; |
| 428 | int outputFormat = 2; |
| 429 | int substitution = 3; |
| 430 | int dualize = 3; |
| 431 | int preSolve = 5; |
| 432 | int doSprint = -1; |
| 433 | int testOsiParameters = -1; |
| 434 | int createSolver = 0; |
| 435 | ClpSimplex * lpSolver; |
| 436 | OsiClpSolverInterface * clpSolver; |
| 437 | if (model_.solver()) { |
| 438 | clpSolver = dynamic_cast<OsiClpSolverInterface *> (model_.solver()); |
| 439 | assert (clpSolver); |
| 440 | lpSolver = clpSolver->getModelPtr(); |
| 441 | assert (lpSolver); |
| 442 | } else { |
| 443 | lpSolver = new ClpSimplex(); |
| 444 | clpSolver = new OsiClpSolverInterface(lpSolver, true); |
| 445 | createSolver = 1 ; |
| 446 | } |
| 447 | parameters_[whichParam(BASISIN, numberParameters_, parameters_)].setStringValue(importBasisFile); |
| 448 | parameters_[whichParam(PRIORITYIN, numberParameters_, parameters_)].setStringValue(importPriorityFile); |
| 449 | parameters_[whichParam(BASISOUT, numberParameters_, parameters_)].setStringValue(exportBasisFile); |
| 450 | parameters_[whichParam(DEBUG, numberParameters_, parameters_)].setStringValue(debugFile); |
| 451 | parameters_[whichParam(PRINTMASK, numberParameters_, parameters_)].setStringValue(printMask); |
| 452 | parameters_[whichParam(DIRECTORY, numberParameters_, parameters_)].setStringValue(directory); |
| 453 | parameters_[whichParam(DIRSAMPLE, numberParameters_, parameters_)].setStringValue(dirSample); |
| 454 | parameters_[whichParam(DIRNETLIB, numberParameters_, parameters_)].setStringValue(dirNetlib); |
| 455 | parameters_[whichParam(DIRMIPLIB, numberParameters_, parameters_)].setStringValue(dirMiplib); |
| 456 | parameters_[whichParam(DUALBOUND, numberParameters_, parameters_)].setDoubleValue(lpSolver->dualBound()); |
| 457 | parameters_[whichParam(DUALTOLERANCE, numberParameters_, parameters_)].setDoubleValue(lpSolver->dualTolerance()); |
| 458 | parameters_[whichParam(EXPORT, numberParameters_, parameters_)].setStringValue(exportFile); |
| 459 | parameters_[whichParam(IDIOT, numberParameters_, parameters_)].setIntValue(doIdiot); |
| 460 | parameters_[whichParam(IMPORT, numberParameters_, parameters_)].setStringValue(importFile); |
| 461 | parameters_[whichParam(PRESOLVETOLERANCE, numberParameters_, parameters_)].setDoubleValue(1.0e-8); |
| 462 | int iParam = whichParam(SOLVERLOGLEVEL, numberParameters_, parameters_); |
| 463 | int value = 1; |
| 464 | clpSolver->messageHandler()->setLogLevel(1) ; |
| 465 | lpSolver->setLogLevel(1); |
| 466 | parameters_[iParam].setIntValue(value); |
| 467 | iParam = whichParam(LOGLEVEL, numberParameters_, parameters_); |
| 468 | model_.messageHandler()->setLogLevel(value); |
| 469 | parameters_[iParam].setIntValue(value); |
| 470 | parameters_[whichParam(MAXFACTOR, numberParameters_, parameters_)].setIntValue(lpSolver->factorizationFrequency()); |
| 471 | parameters_[whichParam(MAXITERATION, numberParameters_, parameters_)].setIntValue(lpSolver->maximumIterations()); |
| 472 | parameters_[whichParam(OUTPUTFORMAT, numberParameters_, parameters_)].setIntValue(outputFormat); |
| 473 | parameters_[whichParam(PRESOLVEPASS, numberParameters_, parameters_)].setIntValue(preSolve); |
| 474 | parameters_[whichParam(PERTVALUE, numberParameters_, parameters_)].setIntValue(lpSolver->perturbation()); |
| 475 | parameters_[whichParam(PRIMALTOLERANCE, numberParameters_, parameters_)].setDoubleValue(lpSolver->primalTolerance()); |
| 476 | parameters_[whichParam(PRIMALWEIGHT, numberParameters_, parameters_)].setDoubleValue(lpSolver->infeasibilityCost()); |
| 477 | parameters_[whichParam(RESTORE, numberParameters_, parameters_)].setStringValue(restoreFile); |
| 478 | parameters_[whichParam(SAVE, numberParameters_, parameters_)].setStringValue(saveFile); |
| 479 | //parameters_[whichParam(TIMELIMIT,numberParameters_,parameters_)].setDoubleValue(1.0e8); |
| 480 | parameters_[whichParam(TIMELIMIT_BAB, numberParameters_, parameters_)].setDoubleValue(1.0e8); |
| 481 | parameters_[whichParam(SOLUTION, numberParameters_, parameters_)].setStringValue(solutionFile); |
| 482 | parameters_[whichParam(SAVESOL, numberParameters_, parameters_)].setStringValue(solutionSaveFile); |
| 483 | parameters_[whichParam(SPRINT, numberParameters_, parameters_)].setIntValue(doSprint); |
| 484 | parameters_[whichParam(SUBSTITUTION, numberParameters_, parameters_)].setIntValue(substitution); |
| 485 | parameters_[whichParam(DUALIZE, numberParameters_, parameters_)].setIntValue(dualize); |
| 486 | parameters_[whichParam(NUMBERBEFORE, numberParameters_, parameters_)].setIntValue(model_.numberBeforeTrust()); |
| 487 | parameters_[whichParam(MAXNODES, numberParameters_, parameters_)].setIntValue(model_.getMaximumNodes()); |
| 488 | parameters_[whichParam(STRONGBRANCHING, numberParameters_, parameters_)].setIntValue(model_.numberStrong()); |
| 489 | parameters_[whichParam(INFEASIBILITYWEIGHT, numberParameters_, parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcInfeasibilityWeight)); |
| 490 | parameters_[whichParam(INTEGERTOLERANCE, numberParameters_, parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcIntegerTolerance)); |
| 491 | parameters_[whichParam(INCREMENT, numberParameters_, parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcCutoffIncrement)); |
| 492 | parameters_[whichParam(TESTOSI, numberParameters_, parameters_)].setIntValue(testOsiParameters); |
| 493 | parameters_[whichParam(FPUMPTUNE, numberParameters_, parameters_)].setIntValue(1003); |
| 494 | initialPumpTune = 1003; |
511 | | parameters_[whichParam(ZEROHALFCUTS,numberParameters_,parameters_)].setCurrentOption("off"); |
512 | | #endif |
513 | | parameters_[whichParam(REDSPLITCUTS,numberParameters_,parameters_)].setCurrentOption("off"); |
514 | | parameters_[whichParam(CLIQUECUTS,numberParameters_,parameters_)].setCurrentOption("ifmove"); |
515 | | parameters_[whichParam(MIXEDCUTS,numberParameters_,parameters_)].setCurrentOption("ifmove"); |
516 | | parameters_[whichParam(FLOWCUTS,numberParameters_,parameters_)].setCurrentOption("ifmove"); |
517 | | parameters_[whichParam(TWOMIRCUTS,numberParameters_,parameters_)].setCurrentOption("root"); |
518 | | parameters_[whichParam(LANDPCUTS,numberParameters_,parameters_)].setCurrentOption("off"); |
519 | | parameters_[whichParam(RESIDCUTS,numberParameters_,parameters_)].setCurrentOption("off"); |
520 | | parameters_[whichParam(ROUNDING,numberParameters_,parameters_)].setCurrentOption("on"); |
521 | | parameters_[whichParam(FPUMP,numberParameters_,parameters_)].setCurrentOption("on"); |
522 | | parameters_[whichParam(GREEDY,numberParameters_,parameters_)].setCurrentOption("on"); |
523 | | parameters_[whichParam(COMBINE,numberParameters_,parameters_)].setCurrentOption("on"); |
524 | | parameters_[whichParam(CROSSOVER2,numberParameters_,parameters_)].setCurrentOption("off"); |
525 | | parameters_[whichParam(PIVOTANDCOMPLEMENT,numberParameters_,parameters_)].setCurrentOption("off"); |
526 | | parameters_[whichParam(PIVOTANDFIX,numberParameters_,parameters_)].setCurrentOption("off"); |
527 | | parameters_[whichParam(RANDROUND,numberParameters_,parameters_)].setCurrentOption("off"); |
528 | | parameters_[whichParam(NAIVE,numberParameters_,parameters_)].setCurrentOption("off"); |
529 | | parameters_[whichParam(RINS,numberParameters_,parameters_)].setCurrentOption("off"); |
530 | | parameters_[whichParam(DINS,numberParameters_,parameters_)].setCurrentOption("off"); |
531 | | parameters_[whichParam(RENS,numberParameters_,parameters_)].setCurrentOption("off"); |
532 | | parameters_[whichParam(LOCALTREE,numberParameters_,parameters_)].setCurrentOption("off"); |
533 | | parameters_[whichParam(COSTSTRATEGY,numberParameters_,parameters_)].setCurrentOption("off"); |
534 | | if (createSolver) |
535 | | delete clpSolver; |
| 511 | parameters_[whichParam(ZEROHALFCUTS, numberParameters_, parameters_)].setCurrentOption("off"); |
| 512 | #endif |
| 513 | parameters_[whichParam(REDSPLITCUTS, numberParameters_, parameters_)].setCurrentOption("off"); |
| 514 | parameters_[whichParam(CLIQUECUTS, numberParameters_, parameters_)].setCurrentOption("ifmove"); |
| 515 | parameters_[whichParam(MIXEDCUTS, numberParameters_, parameters_)].setCurrentOption("ifmove"); |
| 516 | parameters_[whichParam(FLOWCUTS, numberParameters_, parameters_)].setCurrentOption("ifmove"); |
| 517 | parameters_[whichParam(TWOMIRCUTS, numberParameters_, parameters_)].setCurrentOption("root"); |
| 518 | parameters_[whichParam(LANDPCUTS, numberParameters_, parameters_)].setCurrentOption("off"); |
| 519 | parameters_[whichParam(RESIDCUTS, numberParameters_, parameters_)].setCurrentOption("off"); |
| 520 | parameters_[whichParam(ROUNDING, numberParameters_, parameters_)].setCurrentOption("on"); |
| 521 | parameters_[whichParam(FPUMP, numberParameters_, parameters_)].setCurrentOption("on"); |
| 522 | parameters_[whichParam(GREEDY, numberParameters_, parameters_)].setCurrentOption("on"); |
| 523 | parameters_[whichParam(COMBINE, numberParameters_, parameters_)].setCurrentOption("on"); |
| 524 | parameters_[whichParam(CROSSOVER2, numberParameters_, parameters_)].setCurrentOption("off"); |
| 525 | parameters_[whichParam(PIVOTANDCOMPLEMENT, numberParameters_, parameters_)].setCurrentOption("off"); |
| 526 | parameters_[whichParam(PIVOTANDFIX, numberParameters_, parameters_)].setCurrentOption("off"); |
| 527 | parameters_[whichParam(RANDROUND, numberParameters_, parameters_)].setCurrentOption("off"); |
| 528 | parameters_[whichParam(NAIVE, numberParameters_, parameters_)].setCurrentOption("off"); |
| 529 | parameters_[whichParam(RINS, numberParameters_, parameters_)].setCurrentOption("off"); |
| 530 | parameters_[whichParam(DINS, numberParameters_, parameters_)].setCurrentOption("off"); |
| 531 | parameters_[whichParam(RENS, numberParameters_, parameters_)].setCurrentOption("off"); |
| 532 | parameters_[whichParam(LOCALTREE, numberParameters_, parameters_)].setCurrentOption("off"); |
| 533 | parameters_[whichParam(COSTSTRATEGY, numberParameters_, parameters_)].setCurrentOption("off"); |
| 534 | if (createSolver) |
| 535 | delete clpSolver; |
539 | | OsiSolverInterface * solver = model_.solver(); |
540 | | OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver); |
541 | | assert (clpSolver); |
542 | | noPrinting_ = (clpSolver->getModelPtr()->logLevel()==0); |
543 | | CoinMessageHandler * generalMessageHandler = clpSolver->messageHandler(); |
544 | | generalMessageHandler->setPrefix(true); |
545 | | ClpSimplex * lpSolver = clpSolver->getModelPtr(); |
546 | | lpSolver->setPerturbation(50); |
547 | | lpSolver->messageHandler()->setPrefix(false); |
548 | | parameters_[whichParam(DUALBOUND,numberParameters_,parameters_)].setDoubleValue(lpSolver->dualBound()); |
549 | | parameters_[whichParam(DUALTOLERANCE,numberParameters_,parameters_)].setDoubleValue(lpSolver->dualTolerance()); |
550 | | int iParam = whichParam(SOLVERLOGLEVEL,numberParameters_,parameters_); |
551 | | int value=parameters_[iParam].intValue(); |
552 | | clpSolver->messageHandler()->setLogLevel(value) ; |
553 | | lpSolver->setLogLevel(value); |
554 | | iParam = whichParam(LOGLEVEL,numberParameters_,parameters_); |
555 | | value=parameters_[iParam].intValue(); |
556 | | model_.messageHandler()->setLogLevel(value); |
557 | | parameters_[whichParam(LOGLEVEL,numberParameters_,parameters_)].setIntValue(model_.logLevel()); |
558 | | parameters_[whichParam(SOLVERLOGLEVEL,numberParameters_,parameters_)].setIntValue(lpSolver->logLevel()); |
559 | | parameters_[whichParam(MAXFACTOR,numberParameters_,parameters_)].setIntValue(lpSolver->factorizationFrequency()); |
560 | | parameters_[whichParam(MAXITERATION,numberParameters_,parameters_)].setIntValue(lpSolver->maximumIterations()); |
561 | | parameters_[whichParam(PERTVALUE,numberParameters_,parameters_)].setIntValue(lpSolver->perturbation()); |
562 | | parameters_[whichParam(PRIMALTOLERANCE,numberParameters_,parameters_)].setDoubleValue(lpSolver->primalTolerance()); |
563 | | parameters_[whichParam(PRIMALWEIGHT,numberParameters_,parameters_)].setDoubleValue(lpSolver->infeasibilityCost()); |
564 | | parameters_[whichParam(NUMBERBEFORE,numberParameters_,parameters_)].setIntValue(model_.numberBeforeTrust()); |
565 | | parameters_[whichParam(MAXNODES,numberParameters_,parameters_)].setIntValue(model_.getMaximumNodes()); |
566 | | parameters_[whichParam(STRONGBRANCHING,numberParameters_,parameters_)].setIntValue(model_.numberStrong()); |
567 | | parameters_[whichParam(INFEASIBILITYWEIGHT,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcInfeasibilityWeight)); |
568 | | parameters_[whichParam(INTEGERTOLERANCE,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcIntegerTolerance)); |
569 | | parameters_[whichParam(INCREMENT,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcCutoffIncrement)); |
| 539 | OsiSolverInterface * solver = model_.solver(); |
| 540 | OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver); |
| 541 | assert (clpSolver); |
| 542 | noPrinting_ = (clpSolver->getModelPtr()->logLevel() == 0); |
| 543 | CoinMessageHandler * generalMessageHandler = clpSolver->messageHandler(); |
| 544 | generalMessageHandler->setPrefix(true); |
| 545 | ClpSimplex * lpSolver = clpSolver->getModelPtr(); |
| 546 | lpSolver->setPerturbation(50); |
| 547 | lpSolver->messageHandler()->setPrefix(false); |
| 548 | parameters_[whichParam(DUALBOUND, numberParameters_, parameters_)].setDoubleValue(lpSolver->dualBound()); |
| 549 | parameters_[whichParam(DUALTOLERANCE, numberParameters_, parameters_)].setDoubleValue(lpSolver->dualTolerance()); |
| 550 | int iParam = whichParam(SOLVERLOGLEVEL, numberParameters_, parameters_); |
| 551 | int value = parameters_[iParam].intValue(); |
| 552 | clpSolver->messageHandler()->setLogLevel(value) ; |
| 553 | lpSolver->setLogLevel(value); |
| 554 | iParam = whichParam(LOGLEVEL, numberParameters_, parameters_); |
| 555 | value = parameters_[iParam].intValue(); |
| 556 | model_.messageHandler()->setLogLevel(value); |
| 557 | parameters_[whichParam(LOGLEVEL, numberParameters_, parameters_)].setIntValue(model_.logLevel()); |
| 558 | parameters_[whichParam(SOLVERLOGLEVEL, numberParameters_, parameters_)].setIntValue(lpSolver->logLevel()); |
| 559 | parameters_[whichParam(MAXFACTOR, numberParameters_, parameters_)].setIntValue(lpSolver->factorizationFrequency()); |
| 560 | parameters_[whichParam(MAXITERATION, numberParameters_, parameters_)].setIntValue(lpSolver->maximumIterations()); |
| 561 | parameters_[whichParam(PERTVALUE, numberParameters_, parameters_)].setIntValue(lpSolver->perturbation()); |
| 562 | parameters_[whichParam(PRIMALTOLERANCE, numberParameters_, parameters_)].setDoubleValue(lpSolver->primalTolerance()); |
| 563 | parameters_[whichParam(PRIMALWEIGHT, numberParameters_, parameters_)].setDoubleValue(lpSolver->infeasibilityCost()); |
| 564 | parameters_[whichParam(NUMBERBEFORE, numberParameters_, parameters_)].setIntValue(model_.numberBeforeTrust()); |
| 565 | parameters_[whichParam(MAXNODES, numberParameters_, parameters_)].setIntValue(model_.getMaximumNodes()); |
| 566 | parameters_[whichParam(STRONGBRANCHING, numberParameters_, parameters_)].setIntValue(model_.numberStrong()); |
| 567 | parameters_[whichParam(INFEASIBILITYWEIGHT, numberParameters_, parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcInfeasibilityWeight)); |
| 568 | parameters_[whichParam(INTEGERTOLERANCE, numberParameters_, parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcIntegerTolerance)); |
| 569 | parameters_[whichParam(INCREMENT, numberParameters_, parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcCutoffIncrement)); |
662 | | if (!returnMode) |
663 | | return; |
664 | | if (returnMode==2&&babModel_) |
665 | | model_ = *babModel_; |
666 | | if (model2) { |
667 | | // Only continuous valid |
668 | | // update with basis etc |
669 | | // map states |
670 | | /* clp status |
671 | | -1 - unknown e.g. before solve or if postSolve says not optimal |
672 | | 0 - optimal |
673 | | 1 - primal infeasible |
674 | | 2 - dual infeasible |
675 | | 3 - stopped on iterations or time |
676 | | 4 - stopped due to errors |
677 | | 5 - stopped by event handler (virtual int ClpEventHandler::event()) */ |
678 | | /* cbc status |
679 | | -1 before branchAndBound |
680 | | 0 finished - check isProvenOptimal or isProvenInfeasible to see if solution found |
681 | | (or check value of best solution) |
682 | | 1 stopped - on maxnodes, maxsols, maxtime |
683 | | 2 difficulties so run was abandoned |
684 | | (5 event user programmed event occurred) */ |
685 | | /* clp secondary status of problem - may get extended |
686 | | 0 - none |
687 | | 1 - primal infeasible because dual limit reached OR probably primal |
688 | | infeasible but can't prove it (main status 4) |
689 | | 2 - scaled problem optimal - unscaled problem has primal infeasibilities |
690 | | 3 - scaled problem optimal - unscaled problem has dual infeasibilities |
691 | | 4 - scaled problem optimal - unscaled problem has primal and dual infeasibilities |
692 | | 5 - giving up in primal with flagged variables |
693 | | 6 - failed due to empty problem check |
694 | | 7 - postSolve says not optimal |
695 | | 8 - failed due to bad element check |
696 | | 9 - status was 3 and stopped on time |
697 | | 100 up - translation of enum from ClpEventHandler |
698 | | */ |
699 | | /* cbc secondary status of problem |
700 | | -1 unset (status_ will also be -1) |
701 | | 0 search completed with solution |
702 | | 1 linear relaxation not feasible (or worse than cutoff) |
703 | | 2 stopped on gap |
704 | | 3 stopped on nodes |
705 | | 4 stopped on time |
706 | | 5 stopped on user event |
707 | | 6 stopped on solutions |
708 | | 7 linear relaxation unbounded |
709 | | */ |
710 | | int iStatus = model2->status(); |
711 | | int iStatus2 = model2->secondaryStatus(); |
712 | | if (iStatus==0) { |
713 | | iStatus2=0; |
714 | | } else if (iStatus==1) { |
715 | | iStatus=0; |
716 | | iStatus2=1; // say infeasible |
717 | | } else if (iStatus==2) { |
718 | | iStatus=0; |
719 | | iStatus2=7; // say unbounded |
720 | | } else if (iStatus==3) { |
721 | | iStatus=1; |
722 | | if (iStatus2==9) |
723 | | iStatus2=4; |
724 | | else |
725 | | iStatus2=3; // Use nodes - as closer than solutions |
726 | | } else if (iStatus==4) { |
727 | | iStatus=2; // difficulties |
728 | | iStatus2=0; |
| 662 | if (!returnMode) |
| 663 | return; |
| 664 | if (returnMode == 2 && babModel_) |
| 665 | model_ = *babModel_; |
| 666 | if (model2) { |
| 667 | // Only continuous valid |
| 668 | // update with basis etc |
| 669 | // map states |
| 670 | /* clp status |
| 671 | -1 - unknown e.g. before solve or if postSolve says not optimal |
| 672 | 0 - optimal |
| 673 | 1 - primal infeasible |
| 674 | 2 - dual infeasible |
| 675 | 3 - stopped on iterations or time |
| 676 | 4 - stopped due to errors |
| 677 | 5 - stopped by event handler (virtual int ClpEventHandler::event()) */ |
| 678 | /* cbc status |
| 679 | -1 before branchAndBound |
| 680 | 0 finished - check isProvenOptimal or isProvenInfeasible to see if solution found |
| 681 | (or check value of best solution) |
| 682 | 1 stopped - on maxnodes, maxsols, maxtime |
| 683 | 2 difficulties so run was abandoned |
| 684 | (5 event user programmed event occurred) */ |
| 685 | /* clp secondary status of problem - may get extended |
| 686 | 0 - none |
| 687 | 1 - primal infeasible because dual limit reached OR probably primal |
| 688 | infeasible but can't prove it (main status 4) |
| 689 | 2 - scaled problem optimal - unscaled problem has primal infeasibilities |
| 690 | 3 - scaled problem optimal - unscaled problem has dual infeasibilities |
| 691 | 4 - scaled problem optimal - unscaled problem has primal and dual infeasibilities |
| 692 | 5 - giving up in primal with flagged variables |
| 693 | 6 - failed due to empty problem check |
| 694 | 7 - postSolve says not optimal |
| 695 | 8 - failed due to bad element check |
| 696 | 9 - status was 3 and stopped on time |
| 697 | 100 up - translation of enum from ClpEventHandler |
| 698 | */ |
| 699 | /* cbc secondary status of problem |
| 700 | -1 unset (status_ will also be -1) |
| 701 | 0 search completed with solution |
| 702 | 1 linear relaxation not feasible (or worse than cutoff) |
| 703 | 2 stopped on gap |
| 704 | 3 stopped on nodes |
| 705 | 4 stopped on time |
| 706 | 5 stopped on user event |
| 707 | 6 stopped on solutions |
| 708 | 7 linear relaxation unbounded |
| 709 | */ |
| 710 | int iStatus = model2->status(); |
| 711 | int iStatus2 = model2->secondaryStatus(); |
| 712 | if (iStatus == 0) { |
| 713 | iStatus2 = 0; |
| 714 | } else if (iStatus == 1) { |
| 715 | iStatus = 0; |
| 716 | iStatus2 = 1; // say infeasible |
| 717 | } else if (iStatus == 2) { |
| 718 | iStatus = 0; |
| 719 | iStatus2 = 7; // say unbounded |
| 720 | } else if (iStatus == 3) { |
| 721 | iStatus = 1; |
| 722 | if (iStatus2 == 9) |
| 723 | iStatus2 = 4; |
| 724 | else |
| 725 | iStatus2 = 3; // Use nodes - as closer than solutions |
| 726 | } else if (iStatus == 4) { |
| 727 | iStatus = 2; // difficulties |
| 728 | iStatus2 = 0; |
| 729 | } |
| 730 | model_.setProblemStatus(iStatus); |
| 731 | model_.setSecondaryStatus(iStatus2); |
| 732 | OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (model_.solver()); |
| 733 | ClpSimplex * lpSolver = clpSolver->getModelPtr(); |
| 734 | if (model2 != lpSolver) { |
| 735 | lpSolver->moveInfo(*model2); |
| 736 | } |
| 737 | clpSolver->setWarmStart(NULL); // synchronize bases |
| 738 | if (originalSolver_) { |
| 739 | ClpSimplex * lpSolver2 = originalSolver_->getModelPtr(); |
| 740 | assert (model2 != lpSolver2); |
| 741 | lpSolver2->moveInfo(*model2); |
| 742 | originalSolver_->setWarmStart(NULL); // synchronize bases |
| 743 | } |
| 744 | } else if (returnMode == 1) { |
| 745 | OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (model_.solver()); |
| 746 | ClpSimplex * lpSolver = clpSolver->getModelPtr(); |
| 747 | if (babModel_) { |
| 748 | model_.moveInfo(*babModel_); |
| 749 | int numberColumns = babModel_->getNumCols(); |
| 750 | if (babModel_->bestSolution()) |
| 751 | model_.setBestSolution(babModel_->bestSolution(), numberColumns, babModel_->getMinimizationObjValue()); |
| 752 | OsiClpSolverInterface * clpSolver1 = dynamic_cast< OsiClpSolverInterface*> (babModel_->solver()); |
| 753 | ClpSimplex * lpSolver1 = clpSolver1->getModelPtr(); |
| 754 | if (lpSolver1 != lpSolver && model_.bestSolution()) { |
| 755 | lpSolver->moveInfo(*lpSolver1); |
| 756 | } |
| 757 | } |
| 758 | clpSolver->setWarmStart(NULL); // synchronize bases |
944 | | numberContinuous++; |
945 | | } else { |
946 | | if (fabs(value-floor(value+0.5))>1.0e-12) |
947 | | allIntegerCoeff=false; |
948 | | } |
| 951 | double low = rowLower[iRow]; |
| 952 | if (low > -1.0e20) { |
| 953 | low -= sumFixed; |
| 954 | if (fabs(low - floor(low + 0.5)) > 1.0e-12) |
| 955 | allIntegerCoeff = false; |
| 956 | } |
| 957 | double up = rowUpper[iRow]; |
| 958 | if (up < 1.0e20) { |
| 959 | up -= sumFixed; |
| 960 | if (fabs(up - floor(up + 0.5)) > 1.0e-12) |
| 961 | allIntegerCoeff = false; |
| 962 | } |
| 963 | if (!allIntegerCoeff) |
| 964 | continue; // can't do |
| 965 | if (numberContinuous == 1) { |
| 966 | // see if really integer |
| 967 | // This does not allow for complicated cases |
| 968 | if (low == up) { |
| 969 | if (fabs(value1) > 1.0e-3) { |
| 970 | value1 = 1.0 / value1; |
| 971 | if (fabs(value1 - floor(value1 + 0.5)) < 1.0e-12) { |
| 972 | // integer |
| 973 | changed[numberChanged++] = jColumn1; |
| 974 | solver->setInteger(jColumn1); |
| 975 | if (upper[jColumn1] > 1.0e20) |
| 976 | solver->setColUpper(jColumn1, 1.0e20); |
| 977 | if (lower[jColumn1] < -1.0e20) |
| 978 | solver->setColLower(jColumn1, -1.0e20); |
| 979 | } |
| 980 | } |
| 981 | } else { |
| 982 | if (fabs(value1) > 1.0e-3) { |
| 983 | value1 = 1.0 / value1; |
| 984 | if (fabs(value1 - floor(value1 + 0.5)) < 1.0e-12) { |
| 985 | // This constraint will not stop it being integer |
| 986 | ignore[iRow] = 1; |
| 987 | } |
| 988 | } |
| 989 | } |
| 990 | } else if (numberContinuous == 2) { |
| 991 | if (low == up) { |
| 992 | /* need general theory - for now just look at 2 cases - |
| 993 | 1 - +- 1 one in column and just costs i.e. matching objective |
| 994 | 2 - +- 1 two in column but feeds into G/L row which will try and minimize |
| 995 | */ |
| 996 | if (fabs(value1) == 1.0 && value1*value2 == -1.0 && !lower[jColumn1] |
| 997 | && !lower[jColumn2]) { |
| 998 | int n = 0; |
| 999 | int i; |
| 1000 | double objChange = direction * (objective[jColumn1] + objective[jColumn2]); |
| 1001 | double bound = CoinMin(upper[jColumn1], upper[jColumn2]); |
| 1002 | bound = CoinMin(bound, 1.0e20); |
| 1003 | for ( i = columnStart[jColumn1]; i < columnStart[jColumn1] + columnLength[jColumn1]; i++) { |
| 1004 | int jRow = row[i]; |
| 1005 | double value = element[i]; |
| 1006 | if (jRow != iRow) { |
| 1007 | which[n++] = jRow; |
| 1008 | changeRhs[jRow] = value; |
| 1009 | } |
| 1010 | } |
| 1011 | for ( i = columnStart[jColumn1]; i < columnStart[jColumn1] + columnLength[jColumn1]; i++) { |
| 1012 | int jRow = row[i]; |
| 1013 | double value = element[i]; |
| 1014 | if (jRow != iRow) { |
| 1015 | if (!changeRhs[jRow]) { |
| 1016 | which[n++] = jRow; |
| 1017 | changeRhs[jRow] = value; |
| 1018 | } else { |
| 1019 | changeRhs[jRow] += value; |
| 1020 | } |
| 1021 | } |
| 1022 | } |
| 1023 | if (objChange >= 0.0) { |
| 1024 | // see if all rows OK |
| 1025 | bool good = true; |
| 1026 | for (i = 0; i < n; i++) { |
| 1027 | int jRow = which[i]; |
| 1028 | double value = changeRhs[jRow]; |
| 1029 | if (value) { |
| 1030 | value *= bound; |
| 1031 | if (rowLength[jRow] == 1) { |
| 1032 | if (value > 0.0) { |
| 1033 | double rhs = rowLower[jRow]; |
| 1034 | if (rhs > 0.0) { |
| 1035 | double ratio = rhs / value; |
| 1036 | if (fabs(ratio - floor(ratio + 0.5)) > 1.0e-12) |
| 1037 | good = false; |
| 1038 | } |
| 1039 | } else { |
| 1040 | double rhs = rowUpper[jRow]; |
| 1041 | if (rhs < 0.0) { |
| 1042 | double ratio = rhs / value; |
| 1043 | if (fabs(ratio - floor(ratio + 0.5)) > 1.0e-12) |
| 1044 | good = false; |
| 1045 | } |
| 1046 | } |
| 1047 | } else if (rowLength[jRow] == 2) { |
| 1048 | if (value > 0.0) { |
| 1049 | if (rowLower[jRow] > -1.0e20) |
| 1050 | good = false; |
| 1051 | } else { |
| 1052 | if (rowUpper[jRow] < 1.0e20) |
| 1053 | good = false; |
| 1054 | } |
| 1055 | } else { |
| 1056 | good = false; |
| 1057 | } |
| 1058 | } |
| 1059 | } |
| 1060 | if (good) { |
| 1061 | // both can be integer |
| 1062 | changed[numberChanged++] = jColumn1; |
| 1063 | solver->setInteger(jColumn1); |
| 1064 | if (upper[jColumn1] > 1.0e20) |
| 1065 | solver->setColUpper(jColumn1, 1.0e20); |
| 1066 | if (lower[jColumn1] < -1.0e20) |
| 1067 | solver->setColLower(jColumn1, -1.0e20); |
| 1068 | changed[numberChanged++] = jColumn2; |
| 1069 | solver->setInteger(jColumn2); |
| 1070 | if (upper[jColumn2] > 1.0e20) |
| 1071 | solver->setColUpper(jColumn2, 1.0e20); |
| 1072 | if (lower[jColumn2] < -1.0e20) |
| 1073 | solver->setColLower(jColumn2, -1.0e20); |
| 1074 | } |
| 1075 | } |
| 1076 | // clear |
| 1077 | for (i = 0; i < n; i++) { |
| 1078 | changeRhs[which[i]] = 0.0; |
| 1079 | } |
| 1080 | } |
| 1081 | } |
| 1082 | } |
| 1083 | } |
| 1084 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 1085 | if (upper[iColumn] > lower[iColumn] + 1.0e-8 && !solver->isInteger(iColumn)) { |
| 1086 | double value; |
| 1087 | value = upper[iColumn]; |
| 1088 | if (value < 1.0e20 && fabs(value - floor(value + 0.5)) > 1.0e-12) |
| 1089 | continue; |
| 1090 | value = lower[iColumn]; |
| 1091 | if (value > -1.0e20 && fabs(value - floor(value + 0.5)) > 1.0e-12) |
| 1092 | continue; |
| 1093 | bool integer = true; |
| 1094 | for (CoinBigIndex j = columnStart[iColumn]; j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
| 1095 | int iRow = row[j]; |
| 1096 | if (!ignore[iRow]) { |
| 1097 | integer = false; |
| 1098 | break; |
| 1099 | } |
| 1100 | } |
| 1101 | if (integer) { |
| 1102 | // integer |
| 1103 | changed[numberChanged++] = iColumn; |
| 1104 | solver->setInteger(iColumn); |
| 1105 | if (upper[iColumn] > 1.0e20) |
| 1106 | solver->setColUpper(iColumn, 1.0e20); |
| 1107 | if (lower[iColumn] < -1.0e20) |
| 1108 | solver->setColLower(iColumn, -1.0e20); |
| 1109 | } |
| 1110 | } |
| 1111 | } |
| 1112 | finished = numberChanged == saveNumberChanged; |
| 1113 | } |
| 1114 | delete [] which; |
| 1115 | delete [] changeRhs; |
| 1116 | delete [] ignore; |
| 1117 | //if (numberInteger&&!noPrinting_) |
| 1118 | //printf("%d integer variables",numberInteger); |
| 1119 | if (changeInt) { |
| 1120 | //if (!noPrinting_) { |
| 1121 | //if (numberChanged) |
| 1122 | // printf(" and %d variables made integer\n",numberChanged); |
| 1123 | //else |
| 1124 | // printf("\n"); |
| 1125 | //} |
| 1126 | delete [] ignore; |
| 1127 | //increment=0.0; |
| 1128 | if (!numberChanged) { |
| 1129 | delete [] changed; |
| 1130 | delete solver; |
| 1131 | return NULL; |
1238 | | if (doAction==11&&!lastSolution) |
1239 | | lastSolution = model.bestSolution(); |
1240 | | assert (((doAction>=0&&doAction<=3)&&!lastSolution)||(doAction==11&&lastSolution)); |
1241 | | double fractionIntFixed = dextra[3]; |
1242 | | double fractionFixed = dextra[4]; |
1243 | | double fixAbove = dextra[2]; |
1244 | | double fixAboveValue = (dextra[5]>0.0) ? dextra[5] : 1.0; |
1245 | | double time1 = CoinCpuTime(); |
1246 | | int leaveIntFree = extra[1]; |
1247 | | OsiSolverInterface * originalSolver = model.solver(); |
1248 | | OsiClpSolverInterface * originalClpSolver = dynamic_cast< OsiClpSolverInterface*> (originalSolver); |
1249 | | ClpSimplex * originalLpSolver = originalClpSolver->getModelPtr(); |
1250 | | int * originalColumns=NULL; |
1251 | | OsiClpSolverInterface * clpSolver; |
1252 | | ClpSimplex * lpSolver; |
1253 | | ClpPresolve pinfo; |
1254 | | assert(originalSolver->getObjSense()>0); |
1255 | | if (doAction==2||doAction==3) { |
1256 | | double * saveLB=NULL; |
1257 | | double * saveUB=NULL; |
1258 | | int numberColumns = originalLpSolver->numberColumns(); |
1259 | | if (fixAbove>0.0) { |
1260 | | double time1 = CoinCpuTime(); |
1261 | | originalClpSolver->initialSolve(); |
1262 | | printf("first solve took %g seconds\n",CoinCpuTime()-time1); |
1263 | | double * columnLower = originalLpSolver->columnLower() ; |
1264 | | double * columnUpper = originalLpSolver->columnUpper() ; |
1265 | | const double * solution = originalLpSolver->primalColumnSolution(); |
1266 | | saveLB = CoinCopyOfArray(columnLower,numberColumns); |
1267 | | saveUB = CoinCopyOfArray(columnUpper,numberColumns); |
1268 | | const double * objective = originalLpSolver->getObjCoefficients() ; |
1269 | | int iColumn; |
1270 | | int nFix=0; |
1271 | | int nArt=0; |
1272 | | for (iColumn=0;iColumn<numberColumns;iColumn++) { |
1273 | | if (objective[iColumn]>fixAbove) { |
1274 | | if (solution[iColumn]<columnLower[iColumn]+1.0e-8) { |
1275 | | columnUpper[iColumn]=columnLower[iColumn]; |
1276 | | nFix++; |
1277 | | } else { |
1278 | | nArt++; |
1279 | | } |
1280 | | } else if (objective[iColumn]<-fixAbove) { |
1281 | | if (solution[iColumn]>columnUpper[iColumn]-1.0e-8) { |
1282 | | columnLower[iColumn]=columnUpper[iColumn]; |
1283 | | nFix++; |
1284 | | } else { |
1285 | | nArt++; |
1286 | | } |
1287 | | } |
1288 | | } |
1289 | | printf("%d artificials fixed, %d left as in solution\n",nFix,nArt); |
1290 | | lpSolver = pinfo.presolvedModel(*originalLpSolver,1.0e-8,true,10); |
1291 | | if (!lpSolver||doAction==2) { |
1292 | | // take off fixing in original |
1293 | | memcpy(columnLower,saveLB,numberColumns*sizeof(double)); |
1294 | | memcpy(columnUpper,saveUB,numberColumns*sizeof(double)); |
1295 | | } |
1296 | | delete [] saveLB; |
1297 | | delete [] saveUB; |
1298 | | if (!lpSolver) { |
1299 | | // try again |
1300 | | pinfo.destroyPresolve(); |
1301 | | lpSolver = pinfo.presolvedModel(*originalLpSolver,1.0e-8,true,10); |
1302 | | assert (lpSolver); |
1303 | | } |
| 1236 | if (doAction == 11 && !lastSolution) |
| 1237 | lastSolution = model.bestSolution(); |
| 1238 | assert (((doAction >= 0 && doAction <= 3) && !lastSolution) || (doAction == 11 && lastSolution)); |
| 1239 | double fractionIntFixed = dextra[3]; |
| 1240 | double fractionFixed = dextra[4]; |
| 1241 | double fixAbove = dextra[2]; |
| 1242 | double fixAboveValue = (dextra[5] > 0.0) ? dextra[5] : 1.0; |
| 1243 | double time1 = CoinCpuTime(); |
| 1244 | int leaveIntFree = extra[1]; |
| 1245 | OsiSolverInterface * originalSolver = model.solver(); |
| 1246 | OsiClpSolverInterface * originalClpSolver = dynamic_cast< OsiClpSolverInterface*> (originalSolver); |
| 1247 | ClpSimplex * originalLpSolver = originalClpSolver->getModelPtr(); |
| 1248 | int * originalColumns = NULL; |
| 1249 | OsiClpSolverInterface * clpSolver; |
| 1250 | ClpSimplex * lpSolver; |
| 1251 | ClpPresolve pinfo; |
| 1252 | assert(originalSolver->getObjSense() > 0); |
| 1253 | if (doAction == 2 || doAction == 3) { |
| 1254 | double * saveLB = NULL; |
| 1255 | double * saveUB = NULL; |
| 1256 | int numberColumns = originalLpSolver->numberColumns(); |
| 1257 | if (fixAbove > 0.0) { |
| 1258 | double time1 = CoinCpuTime(); |
| 1259 | originalClpSolver->initialSolve(); |
| 1260 | printf("first solve took %g seconds\n", CoinCpuTime() - time1); |
| 1261 | double * columnLower = originalLpSolver->columnLower() ; |
| 1262 | double * columnUpper = originalLpSolver->columnUpper() ; |
| 1263 | const double * solution = originalLpSolver->primalColumnSolution(); |
| 1264 | saveLB = CoinCopyOfArray(columnLower, numberColumns); |
| 1265 | saveUB = CoinCopyOfArray(columnUpper, numberColumns); |
| 1266 | const double * objective = originalLpSolver->getObjCoefficients() ; |
| 1267 | int iColumn; |
| 1268 | int nFix = 0; |
| 1269 | int nArt = 0; |
| 1270 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 1271 | if (objective[iColumn] > fixAbove) { |
| 1272 | if (solution[iColumn] < columnLower[iColumn] + 1.0e-8) { |
| 1273 | columnUpper[iColumn] = columnLower[iColumn]; |
| 1274 | nFix++; |
| 1275 | } else { |
| 1276 | nArt++; |
| 1277 | } |
| 1278 | } else if (objective[iColumn] < -fixAbove) { |
| 1279 | if (solution[iColumn] > columnUpper[iColumn] - 1.0e-8) { |
| 1280 | columnLower[iColumn] = columnUpper[iColumn]; |
| 1281 | nFix++; |
| 1282 | } else { |
| 1283 | nArt++; |
| 1284 | } |
| 1285 | } |
| 1286 | } |
| 1287 | printf("%d artificials fixed, %d left as in solution\n", nFix, nArt); |
| 1288 | lpSolver = pinfo.presolvedModel(*originalLpSolver, 1.0e-8, true, 10); |
| 1289 | if (!lpSolver || doAction == 2) { |
| 1290 | // take off fixing in original |
| 1291 | memcpy(columnLower, saveLB, numberColumns*sizeof(double)); |
| 1292 | memcpy(columnUpper, saveUB, numberColumns*sizeof(double)); |
| 1293 | } |
| 1294 | delete [] saveLB; |
| 1295 | delete [] saveUB; |
| 1296 | if (!lpSolver) { |
| 1297 | // try again |
| 1298 | pinfo.destroyPresolve(); |
| 1299 | lpSolver = pinfo.presolvedModel(*originalLpSolver, 1.0e-8, true, 10); |
| 1300 | assert (lpSolver); |
| 1301 | } |
| 1302 | } else { |
| 1303 | lpSolver = pinfo.presolvedModel(*originalLpSolver, 1.0e-8, true, 10); |
| 1304 | assert (lpSolver); |
| 1305 | } |
| 1306 | clpSolver = new OsiClpSolverInterface(lpSolver, true); |
| 1307 | assert(lpSolver == clpSolver->getModelPtr()); |
| 1308 | numberColumns = lpSolver->numberColumns(); |
| 1309 | originalColumns = CoinCopyOfArray(pinfo.originalColumns(), numberColumns); |
| 1310 | doAction = 1; |
1390 | | if (numberPositive==1&&numberNegative==1) |
1391 | | check[iRow]=-1; // try both |
1392 | | if (numberPositive==1&&rowLower[iRow]>-1.0e20) |
1393 | | check[iRow]=iPositive; |
1394 | | else if (numberNegative==1&&rowUpper[iRow]<1.0e20) |
1395 | | check[iRow]=iNegative; |
1396 | | } |
1397 | | for (iColumn=0;iColumn<numberColumns;iColumn++) { |
1398 | | fix[iColumn]=-1; |
1399 | | if (columnUpper[iColumn] > columnLower[iColumn]+1.0e-8) { |
1400 | | if (clpSolver->isInteger(iColumn)) |
1401 | | numberInteger++; |
1402 | | if (columnLower[iColumn]==0.0) { |
1403 | | bool infeasible=false; |
1404 | | fix[iColumn]=0; |
1405 | | // fake upper bound |
1406 | | double saveUpper = columnUpper[iColumn]; |
1407 | | columnUpper[iColumn]=0.0; |
1408 | | for (CoinBigIndex i=columnStart[iColumn]; |
1409 | | i<columnStart[iColumn]+columnLength[iColumn];i++) { |
1410 | | iRow = row[i]; |
1411 | | if (check[iRow]!=-1&&check[iRow]!=iColumn) |
1412 | | continue; // unlikely |
1413 | | // possible row |
1414 | | int infiniteUpper = 0; |
1415 | | int infiniteLower = 0; |
1416 | | double maximumUp = 0.0; |
1417 | | double maximumDown = 0.0; |
1418 | | double newBound; |
1419 | | CoinBigIndex rStart = rowStart[iRow]; |
1420 | | CoinBigIndex rEnd = rowStart[iRow]+rowLength[iRow]; |
1421 | | CoinBigIndex j; |
1422 | | int kColumn; |
1423 | | // Compute possible lower and upper ranges |
1424 | | for (j = rStart; j < rEnd; ++j) { |
1425 | | double value=elementByRow[j]; |
1426 | | kColumn = column[j]; |
1427 | | if (value > 0.0) { |
1428 | | if (columnUpper[kColumn] >= large) { |
1429 | | ++infiniteUpper; |
1430 | | } else { |
1431 | | maximumUp += columnUpper[kColumn] * value; |
1432 | | } |
1433 | | if (columnLower[kColumn] <= -large) { |
1434 | | ++infiniteLower; |
1435 | | } else { |
1436 | | maximumDown += columnLower[kColumn] * value; |
1437 | | } |
1438 | | } else if (value<0.0) { |
1439 | | if (columnUpper[kColumn] >= large) { |
1440 | | ++infiniteLower; |
1441 | | } else { |
1442 | | maximumDown += columnUpper[kColumn] * value; |
1443 | | } |
1444 | | if (columnLower[kColumn] <= -large) { |
1445 | | ++infiniteUpper; |
1446 | | } else { |
1447 | | maximumUp += columnLower[kColumn] * value; |
1448 | | } |
1449 | | } |
1450 | | } |
1451 | | // Build in a margin of error |
1452 | | maximumUp += 1.0e-8*fabs(maximumUp); |
1453 | | maximumDown -= 1.0e-8*fabs(maximumDown); |
1454 | | double maxUp = maximumUp+infiniteUpper*1.0e31; |
1455 | | double maxDown = maximumDown-infiniteLower*1.0e31; |
1456 | | if (maxUp <= rowUpper[iRow] + tolerance && |
1457 | | maxDown >= rowLower[iRow] - tolerance) { |
1458 | | //printf("Redundant row in vubs %d\n",iRow); |
1459 | | } else { |
1460 | | if (maxUp < rowLower[iRow] -100.0*tolerance || |
1461 | | maxDown > rowUpper[iRow]+100.0*tolerance) { |
1462 | | infeasible=true; |
1463 | | break; |
1464 | | } |
1465 | | double lower = rowLower[iRow]; |
1466 | | double upper = rowUpper[iRow]; |
1467 | | for (j = rStart; j < rEnd; ++j) { |
1468 | | double value=elementByRow[j]; |
1469 | | kColumn = column[j]; |
1470 | | double nowLower = columnLower[kColumn]; |
1471 | | double nowUpper = columnUpper[kColumn]; |
1472 | | if (value > 0.0) { |
1473 | | // positive value |
1474 | | if (lower>-large) { |
1475 | | if (!infiniteUpper) { |
1476 | | assert(nowUpper < large2); |
1477 | | newBound = nowUpper + |
1478 | | (lower - maximumUp) / value; |
1479 | | // relax if original was large |
1480 | | if (fabs(maximumUp)>1.0e8) |
1481 | | newBound -= 1.0e-12*fabs(maximumUp); |
1482 | | } else if (infiniteUpper==1&&nowUpper>large) { |
1483 | | newBound = (lower -maximumUp) / value; |
1484 | | // relax if original was large |
1485 | | if (fabs(maximumUp)>1.0e8) |
1486 | | newBound -= 1.0e-12*fabs(maximumUp); |
1487 | | } else { |
1488 | | newBound = -COIN_DBL_MAX; |
1489 | | } |
1490 | | if (newBound > nowLower + 1.0e-12&&newBound>-large) { |
1491 | | // Tighten the lower bound |
1492 | | // check infeasible (relaxed) |
1493 | | if (nowUpper < newBound) { |
1494 | | if (nowUpper - newBound < |
1495 | | -100.0*tolerance) { |
1496 | | infeasible=true; |
1497 | | break; |
1498 | | } |
1499 | | } |
1500 | | } |
1501 | | } |
1502 | | if (upper <large) { |
1503 | | if (!infiniteLower) { |
1504 | | assert(nowLower >- large2); |
1505 | | newBound = nowLower + |
1506 | | (upper - maximumDown) / value; |
1507 | | // relax if original was large |
1508 | | if (fabs(maximumDown)>1.0e8) |
1509 | | newBound += 1.0e-12*fabs(maximumDown); |
1510 | | } else if (infiniteLower==1&&nowLower<-large) { |
1511 | | newBound = (upper - maximumDown) / value; |
1512 | | // relax if original was large |
1513 | | if (fabs(maximumDown)>1.0e8) |
1514 | | newBound += 1.0e-12*fabs(maximumDown); |
1515 | | } else { |
1516 | | newBound = COIN_DBL_MAX; |
1517 | | } |
1518 | | if (newBound < nowUpper - 1.0e-12&&newBound<large) { |
1519 | | // Tighten the upper bound |
1520 | | // check infeasible (relaxed) |
1521 | | if (nowLower > newBound) { |
1522 | | if (newBound - nowLower < |
1523 | | -100.0*tolerance) { |
1524 | | infeasible=true; |
1525 | | break; |
1526 | | } else { |
1527 | | newBound=nowLower; |
1528 | | } |
1529 | | } |
1530 | | if (!newBound||(clpSolver->isInteger(kColumn)&&newBound<0.999)) { |
1531 | | // fix to zero |
1532 | | if (!mark[kColumn]) { |
1533 | | otherColumn[numberOther++]=kColumn; |
1534 | | mark[kColumn]=1; |
1535 | | if (check[iRow]==-1) |
1536 | | check[iRow]=iColumn; |
1537 | | else |
1538 | | assert(check[iRow]==iColumn); |
1539 | | } |
1540 | | } |
1541 | | } |
1542 | | } |
1543 | | } else { |
1544 | | // negative value |
1545 | | if (lower>-large) { |
1546 | | if (!infiniteUpper) { |
1547 | | assert(nowLower < large2); |
1548 | | newBound = nowLower + |
1549 | | (lower - maximumUp) / value; |
1550 | | // relax if original was large |
1551 | | if (fabs(maximumUp)>1.0e8) |
1552 | | newBound += 1.0e-12*fabs(maximumUp); |
1553 | | } else if (infiniteUpper==1&&nowLower<-large) { |
1554 | | newBound = (lower -maximumUp) / value; |
1555 | | // relax if original was large |
1556 | | if (fabs(maximumUp)>1.0e8) |
1557 | | newBound += 1.0e-12*fabs(maximumUp); |
1558 | | } else { |
1559 | | newBound = COIN_DBL_MAX; |
1560 | | } |
1561 | | if (newBound < nowUpper - 1.0e-12&&newBound<large) { |
1562 | | // Tighten the upper bound |
1563 | | // check infeasible (relaxed) |
1564 | | if (nowLower > newBound) { |
1565 | | if (newBound - nowLower < |
1566 | | -100.0*tolerance) { |
1567 | | infeasible=true; |
1568 | | break; |
1569 | | } else { |
1570 | | newBound=nowLower; |
1571 | | } |
1572 | | } |
1573 | | if (!newBound||(clpSolver->isInteger(kColumn)&&newBound<0.999)) { |
1574 | | // fix to zero |
1575 | | if (!mark[kColumn]) { |
1576 | | otherColumn[numberOther++]=kColumn; |
1577 | | mark[kColumn]=1; |
1578 | | if (check[iRow]==-1) |
1579 | | check[iRow]=iColumn; |
1580 | | else |
1581 | | assert(check[iRow]==iColumn); |
1582 | | } |
1583 | | } |
1584 | | } |
1585 | | } |
1586 | | if (upper <large) { |
1587 | | if (!infiniteLower) { |
1588 | | assert(nowUpper < large2); |
1589 | | newBound = nowUpper + |
1590 | | (upper - maximumDown) / value; |
1591 | | // relax if original was large |
1592 | | if (fabs(maximumDown)>1.0e8) |
1593 | | newBound -= 1.0e-12*fabs(maximumDown); |
1594 | | } else if (infiniteLower==1&&nowUpper>large) { |
1595 | | newBound = (upper - maximumDown) / value; |
1596 | | // relax if original was large |
1597 | | if (fabs(maximumDown)>1.0e8) |
1598 | | newBound -= 1.0e-12*fabs(maximumDown); |
1599 | | } else { |
1600 | | newBound = -COIN_DBL_MAX; |
1601 | | } |
1602 | | if (newBound > nowLower + 1.0e-12&&newBound>-large) { |
1603 | | // Tighten the lower bound |
1604 | | // check infeasible (relaxed) |
1605 | | if (nowUpper < newBound) { |
1606 | | if (nowUpper - newBound < |
1607 | | -100.0*tolerance) { |
1608 | | infeasible=true; |
1609 | | break; |
1610 | | } |
1611 | | } |
1612 | | } |
1613 | | } |
1614 | | } |
1615 | | } |
1616 | | } |
1617 | | } |
1618 | | for (int i=fixColumn[iColumn];i<numberOther;i++) |
1619 | | mark[otherColumn[i]]=0; |
1620 | | // reset bound unless infeasible |
1621 | | if (!infeasible||!clpSolver->isInteger(iColumn)) |
1622 | | columnUpper[iColumn]=saveUpper; |
1623 | | else if (clpSolver->isInteger(iColumn)) |
1624 | | columnLower[iColumn]=1.0; |
1625 | | } |
| 1395 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 1396 | fix[iColumn] = -1; |
| 1397 | if (columnUpper[iColumn] > columnLower[iColumn] + 1.0e-8) { |
| 1398 | if (clpSolver->isInteger(iColumn)) |
| 1399 | numberInteger++; |
| 1400 | if (columnLower[iColumn] == 0.0) { |
| 1401 | bool infeasible = false; |
| 1402 | fix[iColumn] = 0; |
| 1403 | // fake upper bound |
| 1404 | double saveUpper = columnUpper[iColumn]; |
| 1405 | columnUpper[iColumn] = 0.0; |
| 1406 | for (CoinBigIndex i = columnStart[iColumn]; |
| 1407 | i < columnStart[iColumn] + columnLength[iColumn]; i++) { |
| 1408 | iRow = row[i]; |
| 1409 | if (check[iRow] != -1 && check[iRow] != iColumn) |
| 1410 | continue; // unlikely |
| 1411 | // possible row |
| 1412 | int infiniteUpper = 0; |
| 1413 | int infiniteLower = 0; |
| 1414 | double maximumUp = 0.0; |
| 1415 | double maximumDown = 0.0; |
| 1416 | double newBound; |
| 1417 | CoinBigIndex rStart = rowStart[iRow]; |
| 1418 | CoinBigIndex rEnd = rowStart[iRow] + rowLength[iRow]; |
| 1419 | CoinBigIndex j; |
| 1420 | int kColumn; |
| 1421 | // Compute possible lower and upper ranges |
| 1422 | for (j = rStart; j < rEnd; ++j) { |
| 1423 | double value = elementByRow[j]; |
| 1424 | kColumn = column[j]; |
| 1425 | if (value > 0.0) { |
| 1426 | if (columnUpper[kColumn] >= large) { |
| 1427 | ++infiniteUpper; |
| 1428 | } else { |
| 1429 | maximumUp += columnUpper[kColumn] * value; |
| 1430 | } |
| 1431 | if (columnLower[kColumn] <= -large) { |
| 1432 | ++infiniteLower; |
| 1433 | } else { |
| 1434 | maximumDown += columnLower[kColumn] * value; |
| 1435 | } |
| 1436 | } else if (value < 0.0) { |
| 1437 | if (columnUpper[kColumn] >= large) { |
| 1438 | ++infiniteLower; |
| 1439 | } else { |
| 1440 | maximumDown += columnUpper[kColumn] * value; |
| 1441 | } |
| 1442 | if (columnLower[kColumn] <= -large) { |
| 1443 | ++infiniteUpper; |
| 1444 | } else { |
| 1445 | maximumUp += columnLower[kColumn] * value; |
| 1446 | } |
| 1447 | } |
| 1448 | } |
| 1449 | // Build in a margin of error |
| 1450 | maximumUp += 1.0e-8 * fabs(maximumUp); |
| 1451 | maximumDown -= 1.0e-8 * fabs(maximumDown); |
| 1452 | double maxUp = maximumUp + infiniteUpper * 1.0e31; |
| 1453 | double maxDown = maximumDown - infiniteLower * 1.0e31; |
| 1454 | if (maxUp <= rowUpper[iRow] + tolerance && |
| 1455 | maxDown >= rowLower[iRow] - tolerance) { |
| 1456 | //printf("Redundant row in vubs %d\n",iRow); |
| 1457 | } else { |
| 1458 | if (maxUp < rowLower[iRow] - 100.0*tolerance || |
| 1459 | maxDown > rowUpper[iRow] + 100.0*tolerance) { |
| 1460 | infeasible = true; |
| 1461 | break; |
| 1462 | } |
| 1463 | double lower = rowLower[iRow]; |
| 1464 | double upper = rowUpper[iRow]; |
| 1465 | for (j = rStart; j < rEnd; ++j) { |
| 1466 | double value = elementByRow[j]; |
| 1467 | kColumn = column[j]; |
| 1468 | double nowLower = columnLower[kColumn]; |
| 1469 | double nowUpper = columnUpper[kColumn]; |
| 1470 | if (value > 0.0) { |
| 1471 | // positive value |
| 1472 | if (lower > -large) { |
| 1473 | if (!infiniteUpper) { |
| 1474 | assert(nowUpper < large2); |
| 1475 | newBound = nowUpper + |
| 1476 | (lower - maximumUp) / value; |
| 1477 | // relax if original was large |
| 1478 | if (fabs(maximumUp) > 1.0e8) |
| 1479 | newBound -= 1.0e-12 * fabs(maximumUp); |
| 1480 | } else if (infiniteUpper == 1 && nowUpper > large) { |
| 1481 | newBound = (lower - maximumUp) / value; |
| 1482 | // relax if original was large |
| 1483 | if (fabs(maximumUp) > 1.0e8) |
| 1484 | newBound -= 1.0e-12 * fabs(maximumUp); |
| 1485 | } else { |
| 1486 | newBound = -COIN_DBL_MAX; |
| 1487 | } |
| 1488 | if (newBound > nowLower + 1.0e-12 && newBound > -large) { |
| 1489 | // Tighten the lower bound |
| 1490 | // check infeasible (relaxed) |
| 1491 | if (nowUpper < newBound) { |
| 1492 | if (nowUpper - newBound < |
| 1493 | -100.0*tolerance) { |
| 1494 | infeasible = true; |
| 1495 | break; |
| 1496 | } |
| 1497 | } |
| 1498 | } |
| 1499 | } |
| 1500 | if (upper < large) { |
| 1501 | if (!infiniteLower) { |
| 1502 | assert(nowLower > - large2); |
| 1503 | newBound = nowLower + |
| 1504 | (upper - maximumDown) / value; |
| 1505 | // relax if original was large |
| 1506 | if (fabs(maximumDown) > 1.0e8) |
| 1507 | newBound += 1.0e-12 * fabs(maximumDown); |
| 1508 | } else if (infiniteLower == 1 && nowLower < -large) { |
| 1509 | newBound = (upper - maximumDown) / value; |
| 1510 | // relax if original was large |
| 1511 | if (fabs(maximumDown) > 1.0e8) |
| 1512 | newBound += 1.0e-12 * fabs(maximumDown); |
| 1513 | } else { |
| 1514 | newBound = COIN_DBL_MAX; |
| 1515 | } |
| 1516 | if (newBound < nowUpper - 1.0e-12 && newBound < large) { |
| 1517 | // Tighten the upper bound |
| 1518 | // check infeasible (relaxed) |
| 1519 | if (nowLower > newBound) { |
| 1520 | if (newBound - nowLower < |
| 1521 | -100.0*tolerance) { |
| 1522 | infeasible = true; |
| 1523 | break; |
| 1524 | } else { |
| 1525 | newBound = nowLower; |
| 1526 | } |
| 1527 | } |
| 1528 | if (!newBound || (clpSolver->isInteger(kColumn) && newBound < 0.999)) { |
| 1529 | // fix to zero |
| 1530 | if (!mark[kColumn]) { |
| 1531 | otherColumn[numberOther++] = kColumn; |
| 1532 | mark[kColumn] = 1; |
| 1533 | if (check[iRow] == -1) |
| 1534 | check[iRow] = iColumn; |
| 1535 | else |
| 1536 | assert(check[iRow] == iColumn); |
| 1537 | } |
| 1538 | } |
| 1539 | } |
| 1540 | } |
| 1541 | } else { |
| 1542 | // negative value |
| 1543 | if (lower > -large) { |
| 1544 | if (!infiniteUpper) { |
| 1545 | assert(nowLower < large2); |
| 1546 | newBound = nowLower + |
| 1547 | (lower - maximumUp) / value; |
| 1548 | // relax if original was large |
| 1549 | if (fabs(maximumUp) > 1.0e8) |
| 1550 | newBound += 1.0e-12 * fabs(maximumUp); |
| 1551 | } else if (infiniteUpper == 1 && nowLower < -large) { |
| 1552 | newBound = (lower - maximumUp) / value; |
| 1553 | // relax if original was large |
| 1554 | if (fabs(maximumUp) > 1.0e8) |
| 1555 | newBound += 1.0e-12 * fabs(maximumUp); |
| 1556 | } else { |
| 1557 | newBound = COIN_DBL_MAX; |
| 1558 | } |
| 1559 | if (newBound < nowUpper - 1.0e-12 && newBound < large) { |
| 1560 | // Tighten the upper bound |
| 1561 | // check infeasible (relaxed) |
| 1562 | if (nowLower > newBound) { |
| 1563 | if (newBound - nowLower < |
| 1564 | -100.0*tolerance) { |
| 1565 | infeasible = true; |
| 1566 | break; |
| 1567 | } else { |
| 1568 | newBound = nowLower; |
| 1569 | } |
| 1570 | } |
| 1571 | if (!newBound || (clpSolver->isInteger(kColumn) && newBound < 0.999)) { |
| 1572 | // fix to zero |
| 1573 | if (!mark[kColumn]) { |
| 1574 | otherColumn[numberOther++] = kColumn; |
| 1575 | mark[kColumn] = 1; |
| 1576 | if (check[iRow] == -1) |
| 1577 | check[iRow] = iColumn; |
| 1578 | else |
| 1579 | assert(check[iRow] == iColumn); |
| 1580 | } |
| 1581 | } |
| 1582 | } |
| 1583 | } |
| 1584 | if (upper < large) { |
| 1585 | if (!infiniteLower) { |
| 1586 | assert(nowUpper < large2); |
| 1587 | newBound = nowUpper + |
| 1588 | (upper - maximumDown) / value; |
| 1589 | // relax if original was large |
| 1590 | if (fabs(maximumDown) > 1.0e8) |
| 1591 | newBound -= 1.0e-12 * fabs(maximumDown); |
| 1592 | } else if (infiniteLower == 1 && nowUpper > large) { |
| 1593 | newBound = (upper - maximumDown) / value; |
| 1594 | // relax if original was large |
| 1595 | if (fabs(maximumDown) > 1.0e8) |
| 1596 | newBound -= 1.0e-12 * fabs(maximumDown); |
| 1597 | } else { |
| 1598 | newBound = -COIN_DBL_MAX; |
| 1599 | } |
| 1600 | if (newBound > nowLower + 1.0e-12 && newBound > -large) { |
| 1601 | // Tighten the lower bound |
| 1602 | // check infeasible (relaxed) |
| 1603 | if (nowUpper < newBound) { |
| 1604 | if (nowUpper - newBound < |
| 1605 | -100.0*tolerance) { |
| 1606 | infeasible = true; |
| 1607 | break; |
| 1608 | } |
| 1609 | } |
| 1610 | } |
| 1611 | } |
| 1612 | } |
| 1613 | } |
| 1614 | } |
| 1615 | } |
| 1616 | for (int i = fixColumn[iColumn]; i < numberOther; i++) |
| 1617 | mark[otherColumn[i]] = 0; |
| 1618 | // reset bound unless infeasible |
| 1619 | if (!infeasible || !clpSolver->isInteger(iColumn)) |
| 1620 | columnUpper[iColumn] = saveUpper; |
| 1621 | else if (clpSolver->isInteger(iColumn)) |
| 1622 | columnLower[iColumn] = 1.0; |
| 1623 | } |
| 1624 | } |
| 1625 | fixColumn[iColumn+1] = numberOther; |
1764 | | ClpSimplex models[MAXPROB]; |
1765 | | int pass[MAXPROB]; |
1766 | | int kPass=-1; |
1767 | | int kLayer=0; |
1768 | | int skipZero=0; |
1769 | | if (skipZero2==-1) |
1770 | | skipZero2=40; //-1; |
1771 | | /* 0 fixed to 0 by choice |
1772 | | 1 lb of 1 by choice |
1773 | | 2 fixed to 0 by another |
1774 | | 3 as 2 but this go |
1775 | | -1 free |
1776 | | */ |
1777 | | char * state = new char [numberColumns]; |
1778 | | for (iColumn=0;iColumn<numberColumns;iColumn++) |
1779 | | state[iColumn]=-1; |
1780 | | while (true) { |
1781 | | double largest=-0.1; |
1782 | | double smallest=1.1; |
1783 | | int iLargest=-1; |
1784 | | int iSmallest=-1; |
1785 | | int atZero=0; |
1786 | | int atOne=0; |
1787 | | int toZero=0; |
1788 | | int toOne=0; |
1789 | | int numberFree=0; |
1790 | | int numberGreater=0; |
1791 | | columnLower = lpSolver->columnLower(); |
1792 | | columnUpper = lpSolver->columnUpper(); |
1793 | | fullSolution = lpSolver->primalColumnSolution(); |
1794 | | if (doAction==11) { |
1795 | | { |
1796 | | double * columnLower = lpSolver->columnLower(); |
1797 | | double * columnUpper = lpSolver->columnUpper(); |
1798 | | // lpSolver->dual(); |
1799 | | memcpy(columnLower,saveColumnLower,numberColumns*sizeof(double)); |
1800 | | memcpy(columnUpper,saveColumnUpper,numberColumns*sizeof(double)); |
1801 | | // lpSolver->dual(); |
1802 | | int iColumn; |
1803 | | for (iColumn=0;iColumn<numberColumns;iColumn++) { |
1804 | | if (columnUpper[iColumn] > columnLower[iColumn]+1.0e-8) { |
1805 | | if (clpSolver->isInteger(iColumn)) { |
1806 | | double value = lastSolution[iColumn]; |
1807 | | int iValue = static_cast<int> (value+0.5); |
1808 | | assert (fabs(value-static_cast<double> (iValue))<1.0e-3); |
1809 | | assert (iValue>=columnLower[iColumn]&& |
1810 | | iValue<=columnUpper[iColumn]); |
1811 | | columnLower[iColumn]=iValue; |
1812 | | columnUpper[iColumn]=iValue; |
1813 | | } |
1814 | | } |
1815 | | } |
1816 | | lpSolver->initialSolve(solveOptions); |
1817 | | memcpy(columnLower,saveColumnLower,numberColumns*sizeof(double)); |
1818 | | memcpy(columnUpper,saveColumnUpper,numberColumns*sizeof(double)); |
1819 | | } |
1820 | | for (iColumn=0;iColumn<numberColumns;iColumn++) { |
1821 | | if (columnUpper[iColumn] > columnLower[iColumn]+1.0e-8) { |
1822 | | if (clpSolver->isInteger(iColumn)) { |
1823 | | double value = lastSolution[iColumn]; |
1824 | | int iValue = static_cast<int> (value+0.5); |
1825 | | assert (fabs(value-static_cast<double> (iValue))<1.0e-3); |
1826 | | assert (iValue>=columnLower[iColumn]&& |
1827 | | iValue<=columnUpper[iColumn]); |
1828 | | if (!fix[iColumn]) { |
1829 | | if (iValue==0) { |
1830 | | state[iColumn]=0; |
1831 | | assert (!columnLower[iColumn]); |
1832 | | columnUpper[iColumn]=0.0; |
1833 | | } else if (iValue==1) { |
1834 | | state[iColumn]=1; |
1835 | | columnLower[iColumn]=1.0; |
1836 | | } else { |
1837 | | // leave fixed |
1838 | | columnLower[iColumn]=iValue; |
1839 | | columnUpper[iColumn]=iValue; |
1840 | | } |
1841 | | } else if (iValue==0) { |
1842 | | state[iColumn]=10; |
1843 | | columnUpper[iColumn]=0.0; |
1844 | | } else { |
1845 | | // leave fixed |
1846 | | columnLower[iColumn]=iValue; |
1847 | | columnUpper[iColumn]=iValue; |
1848 | | } |
1849 | | } |
1850 | | } |
1851 | | } |
1852 | | int jLayer=0; |
1853 | | int nFixed=-1; |
1854 | | int nTotalFixed=0; |
1855 | | while (nFixed) { |
1856 | | nFixed=0; |
1857 | | for ( iColumn=0;iColumn<numberColumns;iColumn++) { |
1858 | | if (columnUpper[iColumn]==0.0&&fix[iColumn]==jLayer) { |
1859 | | for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) { |
1860 | | int jColumn=otherColumn[i]; |
1861 | | if (columnUpper[jColumn]) { |
1862 | | bool canFix=true; |
1863 | | for (int k=fixColumn2[jColumn];k<fixColumn2[jColumn+1];k++) { |
1864 | | int kColumn=otherColumn2[k]; |
1865 | | if (state[kColumn]==1) { |
1866 | | canFix=false; |
1867 | | break; |
1868 | | } |
1869 | | } |
1870 | | if (canFix) { |
1871 | | columnUpper[jColumn]=0.0; |
1872 | | nFixed++; |
1873 | | } |
1874 | | } |
1875 | | } |
1876 | | } |
1877 | | } |
1878 | | nTotalFixed += nFixed; |
1879 | | jLayer += 100; |
1880 | | } |
1881 | | printf("This fixes %d variables in lower priorities\n",nTotalFixed); |
1882 | | break; |
1883 | | } |
1884 | | for ( iColumn=0;iColumn<numberColumns;iColumn++) { |
1885 | | if (!clpSolver->isInteger(iColumn)||fix[iColumn]>kLayer) |
1886 | | continue; |
1887 | | // skip if fixes nothing |
1888 | | if (fixColumn[iColumn+1]-fixColumn[iColumn]<=skipZero2) |
1889 | | continue; |
1890 | | double value = fullSolution[iColumn]; |
1891 | | if (value>1.00001) { |
1892 | | numberGreater++; |
1893 | | continue; |
1894 | | } |
1895 | | double lower = columnLower[iColumn]; |
1896 | | double upper = columnUpper[iColumn]; |
1897 | | if (lower==upper) { |
1898 | | if (lower) |
1899 | | atOne++; |
1900 | | else |
1901 | | atZero++; |
1902 | | continue; |
1903 | | } |
1904 | | if (value<1.0e-7) { |
1905 | | toZero++; |
1906 | | columnUpper[iColumn]=0.0; |
1907 | | state[iColumn]=10; |
1908 | | continue; |
1909 | | } |
1910 | | if (value>1.0-1.0e-7) { |
1911 | | toOne++; |
1912 | | columnLower[iColumn]=1.0; |
1913 | | state[iColumn]=1; |
1914 | | continue; |
1915 | | } |
1916 | | numberFree++; |
1917 | | // skip if fixes nothing |
1918 | | if (fixColumn[iColumn+1]-fixColumn[iColumn]<=skipZero) |
1919 | | continue; |
1920 | | if (value<smallest) { |
1921 | | smallest=value; |
1922 | | iSmallest=iColumn; |
1923 | | } |
1924 | | if (value>largest) { |
1925 | | largest=value; |
1926 | | iLargest=iColumn; |
1927 | | } |
1928 | | } |
1929 | | if (toZero||toOne) |
1930 | | printf("%d at 0 fixed and %d at one fixed\n",toZero,toOne); |
1931 | | printf("%d variables free, %d fixed to 0, %d to 1 - smallest %g, largest %g\n", |
1932 | | numberFree,atZero,atOne,smallest,largest); |
1933 | | if (numberGreater&&!iPass) |
1934 | | printf("%d variables have value > 1.0\n",numberGreater); |
1935 | | //skipZero2=0; // leave 0 fixing |
1936 | | int jLayer=0; |
1937 | | int nFixed=-1; |
1938 | | int nTotalFixed=0; |
1939 | | while (nFixed) { |
1940 | | nFixed=0; |
1941 | | for ( iColumn=0;iColumn<numberColumns;iColumn++) { |
1942 | | if (columnUpper[iColumn]==0.0&&fix[iColumn]==jLayer) { |
1943 | | for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) { |
1944 | | int jColumn=otherColumn[i]; |
1945 | | if (columnUpper[jColumn]) { |
1946 | | bool canFix=true; |
1947 | | for (int k=fixColumn2[jColumn];k<fixColumn2[jColumn+1];k++) { |
1948 | | int kColumn=otherColumn2[k]; |
1949 | | if (state[kColumn]==1) { |
1950 | | canFix=false; |
1951 | | break; |
1952 | | } |
1953 | | } |
1954 | | if (canFix) { |
1955 | | columnUpper[jColumn]=0.0; |
1956 | | nFixed++; |
1957 | | } |
1958 | | } |
1959 | | } |
1960 | | } |
1961 | | } |
1962 | | nTotalFixed += nFixed; |
1963 | | jLayer += 100; |
1964 | | } |
1965 | | printf("This fixes %d variables in lower priorities\n",nTotalFixed); |
1966 | | if (iLargest<0||numberFree<=leaveIntFree) |
1967 | | break; |
1968 | | double movement; |
1969 | | int way; |
1970 | | if (smallest<=1.0-largest&&smallest<0.2&&largest<fixAboveValue) { |
1971 | | columnUpper[iSmallest]=0.0; |
1972 | | state[iSmallest]=0; |
1973 | | movement=smallest; |
1974 | | way=-1; |
1975 | | } else { |
1976 | | columnLower[iLargest]=1.0; |
1977 | | state[iLargest]=1; |
1978 | | movement=1.0-largest; |
1979 | | way=1; |
1980 | | } |
1981 | | double saveObj = lpSolver->objectiveValue(); |
1982 | | iPass++; |
1983 | | kPass = iPass%MAXPROB; |
1984 | | models[kPass]=*lpSolver; |
1985 | | if (way==-1) { |
1986 | | // fix others |
1987 | | for (int i=fixColumn[iSmallest];i<fixColumn[iSmallest+1];i++) { |
1988 | | int jColumn=otherColumn[i]; |
1989 | | if (state[jColumn]==-1) { |
1990 | | columnUpper[jColumn]=0.0; |
1991 | | state[jColumn]=3; |
1992 | | } |
1993 | | } |
1994 | | } |
1995 | | pass[kPass]=iPass; |
1996 | | double maxCostUp = COIN_DBL_MAX; |
1997 | | objective = lpSolver->getObjCoefficients() ; |
1998 | | if (way==-1) |
1999 | | maxCostUp= (1.0-movement)*objective[iSmallest]; |
2000 | | lpSolver->setDualObjectiveLimit(saveObj+maxCostUp); |
2001 | | crunchIt(lpSolver); |
2002 | | double moveObj = lpSolver->objectiveValue()-saveObj; |
2003 | | printf("movement %s was %g costing %g\n", |
2004 | | (way==-1) ? "down" : "up",movement,moveObj); |
2005 | | if (way==-1&&(moveObj>=maxCostUp||lpSolver->status())) { |
2006 | | // go up |
2007 | | columnLower = models[kPass].columnLower(); |
2008 | | columnUpper = models[kPass].columnUpper(); |
2009 | | columnLower[iSmallest]=1.0; |
2010 | | columnUpper[iSmallest]=saveColumnUpper[iSmallest]; |
2011 | | *lpSolver=models[kPass]; |
2012 | | columnLower = lpSolver->columnLower(); |
2013 | | columnUpper = lpSolver->columnUpper(); |
2014 | | fullSolution = lpSolver->primalColumnSolution(); |
2015 | | dj = lpSolver->dualColumnSolution(); |
2016 | | columnLower[iSmallest]=1.0; |
2017 | | columnUpper[iSmallest]=saveColumnUpper[iSmallest]; |
2018 | | state[iSmallest]=1; |
2019 | | // unfix others |
2020 | | for (int i=fixColumn[iSmallest];i<fixColumn[iSmallest+1];i++) { |
2021 | | int jColumn=otherColumn[i]; |
2022 | | if (state[jColumn]==3) { |
2023 | | columnUpper[jColumn]=saveColumnUpper[jColumn]; |
2024 | | state[jColumn]=-1; |
2025 | | } |
2026 | | } |
2027 | | crunchIt(lpSolver); |
2028 | | } |
2029 | | models[kPass]=*lpSolver; |
| 1762 | ClpSimplex models[MAXPROB]; |
| 1763 | int pass[MAXPROB]; |
| 1764 | int kPass = -1; |
| 1765 | int kLayer = 0; |
| 1766 | int skipZero = 0; |
| 1767 | if (skipZero2 == -1) |
| 1768 | skipZero2 = 40; //-1; |
| 1769 | /* 0 fixed to 0 by choice |
| 1770 | 1 lb of 1 by choice |
| 1771 | 2 fixed to 0 by another |
| 1772 | 3 as 2 but this go |
| 1773 | -1 free |
| 1774 | */ |
| 1775 | char * state = new char [numberColumns]; |
| 1776 | for (iColumn = 0; iColumn < numberColumns; iColumn++) |
| 1777 | state[iColumn] = -1; |
| 1778 | while (true) { |
| 1779 | double largest = -0.1; |
| 1780 | double smallest = 1.1; |
| 1781 | int iLargest = -1; |
| 1782 | int iSmallest = -1; |
| 1783 | int atZero = 0; |
| 1784 | int atOne = 0; |
| 1785 | int toZero = 0; |
| 1786 | int toOne = 0; |
| 1787 | int numberFree = 0; |
| 1788 | int numberGreater = 0; |
| 1789 | columnLower = lpSolver->columnLower(); |
| 1790 | columnUpper = lpSolver->columnUpper(); |
| 1791 | fullSolution = lpSolver->primalColumnSolution(); |
| 1792 | if (doAction == 11) { |
| 1793 | { |
| 1794 | double * columnLower = lpSolver->columnLower(); |
| 1795 | double * columnUpper = lpSolver->columnUpper(); |
| 1796 | // lpSolver->dual(); |
| 1797 | memcpy(columnLower, saveColumnLower, numberColumns*sizeof(double)); |
| 1798 | memcpy(columnUpper, saveColumnUpper, numberColumns*sizeof(double)); |
| 1799 | // lpSolver->dual(); |
| 1800 | int iColumn; |
| 1801 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 1802 | if (columnUpper[iColumn] > columnLower[iColumn] + 1.0e-8) { |
| 1803 | if (clpSolver->isInteger(iColumn)) { |
| 1804 | double value = lastSolution[iColumn]; |
| 1805 | int iValue = static_cast<int> (value + 0.5); |
| 1806 | assert (fabs(value - static_cast<double> (iValue)) < 1.0e-3); |
| 1807 | assert (iValue >= columnLower[iColumn] && |
| 1808 | iValue <= columnUpper[iColumn]); |
| 1809 | columnLower[iColumn] = iValue; |
| 1810 | columnUpper[iColumn] = iValue; |
| 1811 | } |
| 1812 | } |
| 1813 | } |
| 1814 | lpSolver->initialSolve(solveOptions); |
| 1815 | memcpy(columnLower, saveColumnLower, numberColumns*sizeof(double)); |
| 1816 | memcpy(columnUpper, saveColumnUpper, numberColumns*sizeof(double)); |
| 1817 | } |
| 1818 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 1819 | if (columnUpper[iColumn] > columnLower[iColumn] + 1.0e-8) { |
| 1820 | if (clpSolver->isInteger(iColumn)) { |
| 1821 | double value = lastSolution[iColumn]; |
| 1822 | int iValue = static_cast<int> (value + 0.5); |
| 1823 | assert (fabs(value - static_cast<double> (iValue)) < 1.0e-3); |
| 1824 | assert (iValue >= columnLower[iColumn] && |
| 1825 | iValue <= columnUpper[iColumn]); |
| 1826 | if (!fix[iColumn]) { |
| 1827 | if (iValue == 0) { |
| 1828 | state[iColumn] = 0; |
| 1829 | assert (!columnLower[iColumn]); |
| 1830 | columnUpper[iColumn] = 0.0; |
| 1831 | } else if (iValue == 1) { |
| 1832 | state[iColumn] = 1; |
| 1833 | columnLower[iColumn] = 1.0; |
| 1834 | } else { |
| 1835 | // leave fixed |
| 1836 | columnLower[iColumn] = iValue; |
| 1837 | columnUpper[iColumn] = iValue; |
| 1838 | } |
| 1839 | } else if (iValue == 0) { |
| 1840 | state[iColumn] = 10; |
| 1841 | columnUpper[iColumn] = 0.0; |
| 1842 | } else { |
| 1843 | // leave fixed |
| 1844 | columnLower[iColumn] = iValue; |
| 1845 | columnUpper[iColumn] = iValue; |
| 1846 | } |
| 1847 | } |
| 1848 | } |
| 1849 | } |
| 1850 | int jLayer = 0; |
| 1851 | int nFixed = -1; |
| 1852 | int nTotalFixed = 0; |
| 1853 | while (nFixed) { |
| 1854 | nFixed = 0; |
| 1855 | for ( iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 1856 | if (columnUpper[iColumn] == 0.0 && fix[iColumn] == jLayer) { |
| 1857 | for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) { |
| 1858 | int jColumn = otherColumn[i]; |
| 1859 | if (columnUpper[jColumn]) { |
| 1860 | bool canFix = true; |
| 1861 | for (int k = fixColumn2[jColumn]; k < fixColumn2[jColumn+1]; k++) { |
| 1862 | int kColumn = otherColumn2[k]; |
| 1863 | if (state[kColumn] == 1) { |
| 1864 | canFix = false; |
| 1865 | break; |
| 1866 | } |
| 1867 | } |
| 1868 | if (canFix) { |
| 1869 | columnUpper[jColumn] = 0.0; |
| 1870 | nFixed++; |
| 1871 | } |
| 1872 | } |
| 1873 | } |
| 1874 | } |
| 1875 | } |
| 1876 | nTotalFixed += nFixed; |
| 1877 | jLayer += 100; |
| 1878 | } |
| 1879 | printf("This fixes %d variables in lower priorities\n", nTotalFixed); |
| 1880 | break; |
| 1881 | } |
| 1882 | for ( iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 1883 | if (!clpSolver->isInteger(iColumn) || fix[iColumn] > kLayer) |
| 1884 | continue; |
| 1885 | // skip if fixes nothing |
| 1886 | if (fixColumn[iColumn+1] - fixColumn[iColumn] <= skipZero2) |
| 1887 | continue; |
| 1888 | double value = fullSolution[iColumn]; |
| 1889 | if (value > 1.00001) { |
| 1890 | numberGreater++; |
| 1891 | continue; |
| 1892 | } |
| 1893 | double lower = columnLower[iColumn]; |
| 1894 | double upper = columnUpper[iColumn]; |
| 1895 | if (lower == upper) { |
| 1896 | if (lower) |
| 1897 | atOne++; |
| 1898 | else |
| 1899 | atZero++; |
| 1900 | continue; |
| 1901 | } |
| 1902 | if (value < 1.0e-7) { |
| 1903 | toZero++; |
| 1904 | columnUpper[iColumn] = 0.0; |
| 1905 | state[iColumn] = 10; |
| 1906 | continue; |
| 1907 | } |
| 1908 | if (value > 1.0 - 1.0e-7) { |
| 1909 | toOne++; |
| 1910 | columnLower[iColumn] = 1.0; |
| 1911 | state[iColumn] = 1; |
| 1912 | continue; |
| 1913 | } |
| 1914 | numberFree++; |
| 1915 | // skip if fixes nothing |
| 1916 | if (fixColumn[iColumn+1] - fixColumn[iColumn] <= skipZero) |
| 1917 | continue; |
| 1918 | if (value < smallest) { |
| 1919 | smallest = value; |
| 1920 | iSmallest = iColumn; |
| 1921 | } |
| 1922 | if (value > largest) { |
| 1923 | largest = value; |
| 1924 | iLargest = iColumn; |
| 1925 | } |
| 1926 | } |
| 1927 | if (toZero || toOne) |
| 1928 | printf("%d at 0 fixed and %d at one fixed\n", toZero, toOne); |
| 1929 | printf("%d variables free, %d fixed to 0, %d to 1 - smallest %g, largest %g\n", |
| 1930 | numberFree, atZero, atOne, smallest, largest); |
| 1931 | if (numberGreater && !iPass) |
| 1932 | printf("%d variables have value > 1.0\n", numberGreater); |
| 1933 | //skipZero2=0; // leave 0 fixing |
| 1934 | int jLayer = 0; |
| 1935 | int nFixed = -1; |
| 1936 | int nTotalFixed = 0; |
| 1937 | while (nFixed) { |
| 1938 | nFixed = 0; |
| 1939 | for ( iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 1940 | if (columnUpper[iColumn] == 0.0 && fix[iColumn] == jLayer) { |
| 1941 | for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) { |
| 1942 | int jColumn = otherColumn[i]; |
| 1943 | if (columnUpper[jColumn]) { |
| 1944 | bool canFix = true; |
| 1945 | for (int k = fixColumn2[jColumn]; k < fixColumn2[jColumn+1]; k++) { |
| 1946 | int kColumn = otherColumn2[k]; |
| 1947 | if (state[kColumn] == 1) { |
| 1948 | canFix = false; |
| 1949 | break; |
| 1950 | } |
| 1951 | } |
| 1952 | if (canFix) { |
| 1953 | columnUpper[jColumn] = 0.0; |
| 1954 | nFixed++; |
| 1955 | } |
| 1956 | } |
| 1957 | } |
| 1958 | } |
| 1959 | } |
| 1960 | nTotalFixed += nFixed; |
| 1961 | jLayer += 100; |
| 1962 | } |
| 1963 | printf("This fixes %d variables in lower priorities\n", nTotalFixed); |
| 1964 | if (iLargest < 0 || numberFree <= leaveIntFree) |
| 1965 | break; |
| 1966 | double movement; |
| 1967 | int way; |
| 1968 | if (smallest <= 1.0 - largest && smallest < 0.2 && largest < fixAboveValue) { |
| 1969 | columnUpper[iSmallest] = 0.0; |
| 1970 | state[iSmallest] = 0; |
| 1971 | movement = smallest; |
| 1972 | way = -1; |
| 1973 | } else { |
| 1974 | columnLower[iLargest] = 1.0; |
| 1975 | state[iLargest] = 1; |
| 1976 | movement = 1.0 - largest; |
| 1977 | way = 1; |
| 1978 | } |
| 1979 | double saveObj = lpSolver->objectiveValue(); |
| 1980 | iPass++; |
| 1981 | kPass = iPass % MAXPROB; |
| 1982 | models[kPass] = *lpSolver; |
| 1983 | if (way == -1) { |
| 1984 | // fix others |
| 1985 | for (int i = fixColumn[iSmallest]; i < fixColumn[iSmallest+1]; i++) { |
| 1986 | int jColumn = otherColumn[i]; |
| 1987 | if (state[jColumn] == -1) { |
| 1988 | columnUpper[jColumn] = 0.0; |
| 1989 | state[jColumn] = 3; |
| 1990 | } |
| 1991 | } |
| 1992 | } |
| 1993 | pass[kPass] = iPass; |
| 1994 | double maxCostUp = COIN_DBL_MAX; |
| 1995 | objective = lpSolver->getObjCoefficients() ; |
| 1996 | if (way == -1) |
| 1997 | maxCostUp = (1.0 - movement) * objective[iSmallest]; |
| 1998 | lpSolver->setDualObjectiveLimit(saveObj + maxCostUp); |
| 1999 | crunchIt(lpSolver); |
| 2000 | double moveObj = lpSolver->objectiveValue() - saveObj; |
| 2001 | printf("movement %s was %g costing %g\n", |
| 2002 | (way == -1) ? "down" : "up", movement, moveObj); |
| 2003 | if (way == -1 && (moveObj >= maxCostUp || lpSolver->status())) { |
| 2004 | // go up |
| 2005 | columnLower = models[kPass].columnLower(); |
| 2006 | columnUpper = models[kPass].columnUpper(); |
| 2007 | columnLower[iSmallest] = 1.0; |
| 2008 | columnUpper[iSmallest] = saveColumnUpper[iSmallest]; |
| 2009 | *lpSolver = models[kPass]; |
| 2010 | columnLower = lpSolver->columnLower(); |
| 2011 | columnUpper = lpSolver->columnUpper(); |
| 2012 | fullSolution = lpSolver->primalColumnSolution(); |
| 2013 | dj = lpSolver->dualColumnSolution(); |
| 2014 | columnLower[iSmallest] = 1.0; |
| 2015 | columnUpper[iSmallest] = saveColumnUpper[iSmallest]; |
| 2016 | state[iSmallest] = 1; |
| 2017 | // unfix others |
| 2018 | for (int i = fixColumn[iSmallest]; i < fixColumn[iSmallest+1]; i++) { |
| 2019 | int jColumn = otherColumn[i]; |
| 2020 | if (state[jColumn] == 3) { |
| 2021 | columnUpper[jColumn] = saveColumnUpper[jColumn]; |
| 2022 | state[jColumn] = -1; |
| 2023 | } |
| 2024 | } |
| 2025 | crunchIt(lpSolver); |
| 2026 | } |
| 2027 | models[kPass] = *lpSolver; |
| 2028 | } |
| 2029 | lpSolver->dual(); |
| 2030 | printf("Fixing took %g seconds\n", CoinCpuTime() - time1); |
| 2031 | columnLower = lpSolver->columnLower(); |
| 2032 | columnUpper = lpSolver->columnUpper(); |
| 2033 | fullSolution = lpSolver->primalColumnSolution(); |
| 2034 | dj = lpSolver->dualColumnSolution(); |
| 2035 | int * sort = new int[numberColumns]; |
| 2036 | double * dsort = new double[numberColumns]; |
| 2037 | int chunk = 20; |
| 2038 | int iRelax = 0; |
| 2039 | //double fractionFixed=6.0/8.0; |
| 2040 | // relax while lots fixed |
| 2041 | while (true) { |
| 2042 | if (skipZero2 > 10 && doAction < 10) |
| 2043 | break; |
| 2044 | iRelax++; |
| 2045 | int n = 0; |
| 2046 | double sum0 = 0.0; |
| 2047 | double sum00 = 0.0; |
| 2048 | double sum1 = 0.0; |
| 2049 | for ( iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 2050 | if (!clpSolver->isInteger(iColumn) || fix[iColumn] > kLayer) |
| 2051 | continue; |
| 2052 | // skip if fixes nothing |
| 2053 | if (fixColumn[iColumn+1] - fixColumn[iColumn] == 0 && doAction < 10) |
| 2054 | continue; |
| 2055 | double djValue = dj[iColumn]; |
| 2056 | if (state[iColumn] == 1) { |
| 2057 | assert (columnLower[iColumn]); |
| 2058 | assert (fullSolution[iColumn] > 0.1); |
| 2059 | if (djValue > 0.0) { |
| 2060 | //printf("YY dj of %d at %g is %g\n",iColumn,value,djValue); |
| 2061 | sum1 += djValue; |
| 2062 | sort[n] = iColumn; |
| 2063 | dsort[n++] = -djValue; |
| 2064 | } else { |
| 2065 | //printf("dj of %d at %g is %g\n",iColumn,value,djValue); |
| 2066 | } |
| 2067 | } else if (state[iColumn] == 0 || state[iColumn] == 10) { |
| 2068 | assert (fullSolution[iColumn] < 0.1); |
| 2069 | assert (!columnUpper[iColumn]); |
| 2070 | double otherValue = 0.0; |
| 2071 | int nn = 0; |
| 2072 | for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) { |
| 2073 | int jColumn = otherColumn[i]; |
| 2074 | if (columnUpper[jColumn] == 0.0) { |
| 2075 | if (dj[jColumn] < -1.0e-5) { |
| 2076 | nn++; |
| 2077 | otherValue += dj[jColumn]; // really need to look at rest |
| 2078 | } |
| 2079 | } |
| 2080 | } |
| 2081 | if (djValue < -1.0e-2 || otherValue < -1.0e-2) { |
| 2082 | //printf("XX dj of %d at %g is %g - %d out of %d contribute %g\n",iColumn,value,djValue, |
| 2083 | // nn,fixColumn[iColumn+1]-fixColumn[iColumn],otherValue); |
| 2084 | if (djValue < 1.0e-8) { |
| 2085 | sum0 -= djValue; |
| 2086 | sum00 -= otherValue; |
| 2087 | sort[n] = iColumn; |
| 2088 | if (djValue < -1.0e-2) |
| 2089 | dsort[n++] = djValue + otherValue; |
| 2090 | else |
| 2091 | dsort[n++] = djValue + 0.001 * otherValue; |
| 2092 | } |
| 2093 | } else { |
| 2094 | //printf("dj of %d at %g is %g - no contribution from %d\n",iColumn,value,djValue, |
| 2095 | // fixColumn[iColumn+1]-fixColumn[iColumn]); |
| 2096 | } |
| 2097 | } |
| 2098 | } |
| 2099 | CoinSort_2(dsort, dsort + n, sort); |
| 2100 | double * originalColumnLower = saveColumnLower; |
| 2101 | double * originalColumnUpper = saveColumnUpper; |
| 2102 | double * lo = CoinCopyOfArray(columnLower, numberColumns); |
| 2103 | double * up = CoinCopyOfArray(columnUpper, numberColumns); |
| 2104 | for (int k = 0; k < CoinMin(chunk, n); k++) { |
| 2105 | iColumn = sort[k]; |
| 2106 | state[iColumn] = -2; |
| 2107 | } |
| 2108 | memcpy(columnLower, originalColumnLower, numberColumns*sizeof(double)); |
| 2109 | memcpy(columnUpper, originalColumnUpper, numberColumns*sizeof(double)); |
| 2110 | int nFixed = 0; |
| 2111 | int nFixed0 = 0; |
| 2112 | int nFixed1 = 0; |
| 2113 | for ( iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 2114 | if (state[iColumn] == 0 || state[iColumn] == 10) { |
| 2115 | columnUpper[iColumn] = 0.0; |
| 2116 | assert (lo[iColumn] == 0.0); |
| 2117 | nFixed++; |
| 2118 | nFixed0++; |
| 2119 | for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) { |
| 2120 | int jColumn = otherColumn[i]; |
| 2121 | if (columnUpper[jColumn]) { |
| 2122 | bool canFix = true; |
| 2123 | for (int k = fixColumn2[jColumn]; k < fixColumn2[jColumn+1]; k++) { |
| 2124 | int kColumn = otherColumn2[k]; |
| 2125 | if (state[kColumn] == 1 || state[kColumn] == -2) { |
| 2126 | canFix = false; |
| 2127 | break; |
| 2128 | } |
| 2129 | } |
| 2130 | if (canFix) { |
| 2131 | columnUpper[jColumn] = 0.0; |
| 2132 | assert (lo[jColumn] == 0.0); |
| 2133 | nFixed++; |
| 2134 | } |
| 2135 | } |
| 2136 | } |
| 2137 | } else if (state[iColumn] == 1) { |
| 2138 | columnLower[iColumn] = 1.0; |
| 2139 | nFixed1++; |
| 2140 | } |
| 2141 | } |
| 2142 | printf("%d fixed %d orig 0 %d 1\n", nFixed, nFixed0, nFixed1); |
| 2143 | int jLayer = 0; |
| 2144 | nFixed = -1; |
| 2145 | int nTotalFixed = 0; |
| 2146 | while (nFixed) { |
| 2147 | nFixed = 0; |
| 2148 | for ( iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 2149 | if (columnUpper[iColumn] == 0.0 && fix[iColumn] == jLayer) { |
| 2150 | for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) { |
| 2151 | int jColumn = otherColumn[i]; |
| 2152 | if (columnUpper[jColumn]) { |
| 2153 | bool canFix = true; |
| 2154 | for (int k = fixColumn2[jColumn]; k < fixColumn2[jColumn+1]; k++) { |
| 2155 | int kColumn = otherColumn2[k]; |
| 2156 | if (state[kColumn] == 1 || state[kColumn] == -2) { |
| 2157 | canFix = false; |
| 2158 | break; |
| 2159 | } |
| 2160 | } |
| 2161 | if (canFix) { |
| 2162 | columnUpper[jColumn] = 0.0; |
| 2163 | assert (lo[jColumn] == 0.0); |
| 2164 | nFixed++; |
| 2165 | } |
| 2166 | } |
| 2167 | } |
| 2168 | } |
| 2169 | } |
| 2170 | nTotalFixed += nFixed; |
| 2171 | jLayer += 100; |
| 2172 | } |
| 2173 | nFixed = 0; |
| 2174 | int nFixedI = 0; |
| 2175 | for ( iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 2176 | if (columnLower[iColumn] == columnUpper[iColumn]) { |
| 2177 | if (clpSolver->isInteger(iColumn)) |
| 2178 | nFixedI++; |
| 2179 | nFixed++; |
| 2180 | } |
| 2181 | } |
| 2182 | printf("This fixes %d variables in lower priorities - total %d (%d integer) - all target %d, int target %d\n", |
| 2183 | nTotalFixed, nFixed, nFixedI, static_cast<int>(fractionFixed*numberColumns), static_cast<int> (fractionIntFixed*numberInteger)); |
| 2184 | int nBad = 0; |
| 2185 | int nRelax = 0; |
| 2186 | for ( iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 2187 | if (lo[iColumn] < columnLower[iColumn] || |
| 2188 | up[iColumn] > columnUpper[iColumn]) { |
| 2189 | printf("bad %d old %g %g, new %g %g\n", iColumn, lo[iColumn], up[iColumn], |
| 2190 | columnLower[iColumn], columnUpper[iColumn]); |
| 2191 | nBad++; |
| 2192 | } |
| 2193 | if (lo[iColumn] > columnLower[iColumn] || |
| 2194 | up[iColumn] < columnUpper[iColumn]) { |
| 2195 | nRelax++; |
| 2196 | } |
| 2197 | } |
| 2198 | printf("%d relaxed\n", nRelax); |
| 2199 | if (iRelax > 20 && nRelax == chunk) |
| 2200 | nRelax = 0; |
| 2201 | if (iRelax > 50) |
| 2202 | nRelax = 0; |
| 2203 | assert (!nBad); |
| 2204 | delete [] lo; |
| 2205 | delete [] up; |
| 2206 | lpSolver->primal(1); |
| 2207 | if (nFixed < fractionFixed*numberColumns || nFixedI < fractionIntFixed*numberInteger || !nRelax) |
| 2208 | break; |
| 2209 | } |
| 2210 | delete [] state; |
| 2211 | delete [] sort; |
| 2212 | delete [] dsort; |
2031 | | lpSolver->dual(); |
2032 | | printf("Fixing took %g seconds\n",CoinCpuTime()-time1); |
2033 | | columnLower = lpSolver->columnLower(); |
2034 | | columnUpper = lpSolver->columnUpper(); |
2035 | | fullSolution = lpSolver->primalColumnSolution(); |
2036 | | dj = lpSolver->dualColumnSolution(); |
2037 | | int * sort = new int[numberColumns]; |
2038 | | double * dsort = new double[numberColumns]; |
2039 | | int chunk=20; |
2040 | | int iRelax=0; |
2041 | | //double fractionFixed=6.0/8.0; |
2042 | | // relax while lots fixed |
2043 | | while (true) { |
2044 | | if (skipZero2>10&&doAction<10) |
2045 | | break; |
2046 | | iRelax++; |
2047 | | int n=0; |
2048 | | double sum0=0.0; |
2049 | | double sum00=0.0; |
2050 | | double sum1=0.0; |
2051 | | for ( iColumn=0;iColumn<numberColumns;iColumn++) { |
2052 | | if (!clpSolver->isInteger(iColumn)||fix[iColumn]>kLayer) |
2053 | | continue; |
2054 | | // skip if fixes nothing |
2055 | | if (fixColumn[iColumn+1]-fixColumn[iColumn]==0&&doAction<10) |
2056 | | continue; |
2057 | | double djValue = dj[iColumn]; |
2058 | | if (state[iColumn]==1) { |
2059 | | assert (columnLower[iColumn]); |
2060 | | assert (fullSolution[iColumn]>0.1); |
2061 | | if (djValue>0.0) { |
2062 | | //printf("YY dj of %d at %g is %g\n",iColumn,value,djValue); |
2063 | | sum1 += djValue; |
2064 | | sort[n]=iColumn; |
2065 | | dsort[n++]=-djValue; |
2066 | | } else { |
2067 | | //printf("dj of %d at %g is %g\n",iColumn,value,djValue); |
2068 | | } |
2069 | | } else if (state[iColumn]==0||state[iColumn]==10) { |
2070 | | assert (fullSolution[iColumn]<0.1); |
2071 | | assert (!columnUpper[iColumn]); |
2072 | | double otherValue=0.0; |
2073 | | int nn=0; |
2074 | | for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) { |
2075 | | int jColumn=otherColumn[i]; |
2076 | | if (columnUpper[jColumn]==0.0) { |
2077 | | if (dj[jColumn]<-1.0e-5) { |
2078 | | nn++; |
2079 | | otherValue += dj[jColumn]; // really need to look at rest |
2080 | | } |
2081 | | } |
2082 | | } |
2083 | | if (djValue<-1.0e-2||otherValue<-1.0e-2) { |
2084 | | //printf("XX dj of %d at %g is %g - %d out of %d contribute %g\n",iColumn,value,djValue, |
2085 | | // nn,fixColumn[iColumn+1]-fixColumn[iColumn],otherValue); |
2086 | | if (djValue<1.0e-8) { |
2087 | | sum0 -= djValue; |
2088 | | sum00 -= otherValue; |
2089 | | sort[n]=iColumn; |
2090 | | if (djValue<-1.0e-2) |
2091 | | dsort[n++]=djValue+otherValue; |
2092 | | else |
2093 | | dsort[n++]=djValue+0.001*otherValue; |
2094 | | } |
2095 | | } else { |
2096 | | //printf("dj of %d at %g is %g - no contribution from %d\n",iColumn,value,djValue, |
2097 | | // fixColumn[iColumn+1]-fixColumn[iColumn]); |
2098 | | } |
2099 | | } |
2100 | | } |
2101 | | CoinSort_2(dsort,dsort+n,sort); |
2102 | | double * originalColumnLower = saveColumnLower; |
2103 | | double * originalColumnUpper = saveColumnUpper; |
2104 | | double * lo = CoinCopyOfArray(columnLower,numberColumns); |
2105 | | double * up = CoinCopyOfArray(columnUpper,numberColumns); |
2106 | | for (int k=0;k<CoinMin(chunk,n);k++) { |
2107 | | iColumn = sort[k]; |
2108 | | state[iColumn]=-2; |
2109 | | } |
2110 | | memcpy(columnLower,originalColumnLower,numberColumns*sizeof(double)); |
2111 | | memcpy(columnUpper,originalColumnUpper,numberColumns*sizeof(double)); |
2112 | | int nFixed=0; |
2113 | | int nFixed0=0; |
2114 | | int nFixed1=0; |
2115 | | for ( iColumn=0;iColumn<numberColumns;iColumn++) { |
2116 | | if (state[iColumn]==0||state[iColumn]==10) { |
2117 | | columnUpper[iColumn]=0.0; |
2118 | | assert (lo[iColumn]==0.0); |
2119 | | nFixed++; |
2120 | | nFixed0++; |
2121 | | for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) { |
2122 | | int jColumn=otherColumn[i]; |
2123 | | if (columnUpper[jColumn]) { |
2124 | | bool canFix=true; |
2125 | | for (int k=fixColumn2[jColumn];k<fixColumn2[jColumn+1];k++) { |
2126 | | int kColumn=otherColumn2[k]; |
2127 | | if (state[kColumn]==1||state[kColumn]==-2) { |
2128 | | canFix=false; |
2129 | | break; |
2130 | | } |
2131 | | } |
2132 | | if (canFix) { |
2133 | | columnUpper[jColumn]=0.0; |
2134 | | assert (lo[jColumn]==0.0); |
2135 | | nFixed++; |
2136 | | } |
2137 | | } |
2138 | | } |
2139 | | } else if (state[iColumn]==1) { |
2140 | | columnLower[iColumn]=1.0; |
2141 | | nFixed1++; |
2142 | | } |
2143 | | } |
2144 | | printf("%d fixed %d orig 0 %d 1\n",nFixed,nFixed0,nFixed1); |
2145 | | int jLayer=0; |
2146 | | nFixed=-1; |
2147 | | int nTotalFixed=0; |
2148 | | while (nFixed) { |
2149 | | nFixed=0; |
2150 | | for ( iColumn=0;iColumn<numberColumns;iColumn++) { |
2151 | | if (columnUpper[iColumn]==0.0&&fix[iColumn]==jLayer) { |
2152 | | for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) { |
2153 | | int jColumn=otherColumn[i]; |
2154 | | if (columnUpper[jColumn]) { |
2155 | | bool canFix=true; |
2156 | | for (int k=fixColumn2[jColumn];k<fixColumn2[jColumn+1];k++) { |
2157 | | int kColumn=otherColumn2[k]; |
2158 | | if (state[kColumn]==1||state[kColumn]==-2) { |
2159 | | canFix=false; |
2160 | | break; |
2161 | | } |
2162 | | } |
2163 | | if (canFix) { |
2164 | | columnUpper[jColumn]=0.0; |
2165 | | assert (lo[jColumn]==0.0); |
2166 | | nFixed++; |
2167 | | } |
2168 | | } |
2169 | | } |
2170 | | } |
2171 | | } |
2172 | | nTotalFixed += nFixed; |
2173 | | jLayer += 100; |
2174 | | } |
2175 | | nFixed=0; |
2176 | | int nFixedI=0; |
2177 | | for ( iColumn=0;iColumn<numberColumns;iColumn++) { |
2178 | | if (columnLower[iColumn]==columnUpper[iColumn]) { |
2179 | | if (clpSolver->isInteger(iColumn)) |
2180 | | nFixedI++; |
2181 | | nFixed++; |
2182 | | } |
2183 | | } |
2184 | | printf("This fixes %d variables in lower priorities - total %d (%d integer) - all target %d, int target %d\n", |
2185 | | nTotalFixed,nFixed,nFixedI,static_cast<int>(fractionFixed*numberColumns),static_cast<int> (fractionIntFixed*numberInteger)); |
2186 | | int nBad=0; |
2187 | | int nRelax=0; |
2188 | | for ( iColumn=0;iColumn<numberColumns;iColumn++) { |
2189 | | if (lo[iColumn]<columnLower[iColumn]|| |
2190 | | up[iColumn]>columnUpper[iColumn]) { |
2191 | | printf("bad %d old %g %g, new %g %g\n",iColumn,lo[iColumn],up[iColumn], |
2192 | | columnLower[iColumn],columnUpper[iColumn]); |
2193 | | nBad++; |
2194 | | } |
2195 | | if (lo[iColumn]>columnLower[iColumn]|| |
2196 | | up[iColumn]<columnUpper[iColumn]) { |
2197 | | nRelax++; |
2198 | | } |
2199 | | } |
2200 | | printf("%d relaxed\n",nRelax); |
2201 | | if (iRelax>20&&nRelax==chunk) |
2202 | | nRelax=0; |
2203 | | if (iRelax>50) |
2204 | | nRelax=0; |
2205 | | assert (!nBad); |
2206 | | delete [] lo; |
2207 | | delete [] up; |
2208 | | lpSolver->primal(1); |
2209 | | if (nFixed<fractionFixed*numberColumns||nFixedI<fractionIntFixed*numberInteger||!nRelax) |
2210 | | break; |
| 2214 | delete [] fix; |
| 2215 | delete [] fixColumn; |
| 2216 | delete [] otherColumn; |
| 2217 | delete [] otherColumn2; |
| 2218 | delete [] fixColumn2; |
| 2219 | // See if was presolved |
| 2220 | if (originalColumns) { |
| 2221 | columnLower = lpSolver->columnLower(); |
| 2222 | columnUpper = lpSolver->columnUpper(); |
| 2223 | for ( iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 2224 | saveColumnLower[iColumn] = columnLower[iColumn]; |
| 2225 | saveColumnUpper[iColumn] = columnUpper[iColumn]; |
| 2226 | } |
| 2227 | pinfo.postsolve(true); |
| 2228 | columnLower = originalLpSolver->columnLower(); |
| 2229 | columnUpper = originalLpSolver->columnUpper(); |
| 2230 | double * newColumnLower = lpSolver->columnLower(); |
| 2231 | double * newColumnUpper = lpSolver->columnUpper(); |
| 2232 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 2233 | int jColumn = originalColumns[iColumn]; |
| 2234 | columnLower[jColumn] = CoinMax(columnLower[jColumn], newColumnLower[iColumn]); |
| 2235 | columnUpper[jColumn] = CoinMin(columnUpper[jColumn], newColumnUpper[iColumn]); |
| 2236 | } |
| 2237 | numberColumns = originalLpSolver->numberColumns(); |
| 2238 | delete [] originalColumns; |
2315 | | if (SOSPriority<0) |
2316 | | SOSPriority=100000; |
2317 | | } |
2318 | | CoinModel coinModel=*si->coinModel(); |
2319 | | assert(coinModel.numberRows()>0); |
2320 | | tightenedModel = coinModel; |
2321 | | int numberRows = coinModel.numberRows(); |
2322 | | // Mark variables |
2323 | | int * whichKnapsack = new int [numberColumns]; |
2324 | | int iRow,iColumn; |
2325 | | for (iColumn=0;iColumn<numberColumns;iColumn++) |
2326 | | whichKnapsack[iColumn]=-1; |
2327 | | int kRow; |
2328 | | bool badModel=false; |
2329 | | // analyze |
2330 | | if (logLevel>1) { |
2331 | | for (iRow=0;iRow<numberRows;iRow++) { |
2332 | | /* Just obvious one at first |
2333 | | positive non unit coefficients |
2334 | | all integer |
2335 | | positive rowUpper |
2336 | | for now - linear (but further down in code may use nonlinear) |
2337 | | column bounds should be tight |
2338 | | */ |
2339 | | //double lower = coinModel.getRowLower(iRow); |
2340 | | double upper = coinModel.getRowUpper(iRow); |
2341 | | if (upper<1.0e10) { |
2342 | | CoinModelLink triple=coinModel.firstInRow(iRow); |
2343 | | bool possible=true; |
2344 | | int n=0; |
2345 | | int n1=0; |
2346 | | while (triple.column()>=0) { |
2347 | | int iColumn = triple.column(); |
2348 | | const char * el = coinModel.getElementAsString(iRow,iColumn); |
2349 | | if (!strcmp("Numeric",el)) { |
2350 | | if (coinModel.columnLower(iColumn)==coinModel.columnUpper(iColumn)) { |
2351 | | triple=coinModel.next(triple); |
2352 | | continue; // fixed |
2353 | | } |
2354 | | double value=coinModel.getElement(iRow,iColumn); |
2355 | | if (value<0.0) { |
2356 | | possible=false; |
2357 | | } else { |
2358 | | n++; |
2359 | | if (value==1.0) |
2360 | | n1++; |
2361 | | if (coinModel.columnLower(iColumn)<0.0) |
2362 | | possible=false; |
2363 | | if (!coinModel.isInteger(iColumn)) |
2364 | | possible=false; |
2365 | | if (whichKnapsack[iColumn]>=0) |
2366 | | possible=false; |
2367 | | } |
2368 | | } else { |
2369 | | possible=false; // non linear |
2370 | | } |
2371 | | triple=coinModel.next(triple); |
2372 | | } |
2373 | | if (n-n1>1&&possible) { |
2374 | | double lower = coinModel.getRowLower(iRow); |
2375 | | double upper = coinModel.getRowUpper(iRow); |
2376 | | CoinModelLink triple=coinModel.firstInRow(iRow); |
2377 | | while (triple.column()>=0) { |
2378 | | int iColumn = triple.column(); |
2379 | | lower -= coinModel.columnLower(iColumn)*triple.value(); |
2380 | | upper -= coinModel.columnLower(iColumn)*triple.value(); |
2381 | | triple=coinModel.next(triple); |
2382 | | } |
2383 | | printf("%d is possible %g <=",iRow,lower); |
2384 | | // print |
2385 | | triple=coinModel.firstInRow(iRow); |
2386 | | while (triple.column()>=0) { |
2387 | | int iColumn = triple.column(); |
2388 | | if (coinModel.columnLower(iColumn)!=coinModel.columnUpper(iColumn)) |
2389 | | printf(" (%d,el %g up %g)",iColumn,triple.value(), |
2390 | | coinModel.columnUpper(iColumn)-coinModel.columnLower(iColumn)); |
2391 | | triple=coinModel.next(triple); |
2392 | | } |
2393 | | printf(" <= %g\n",upper); |
2394 | | } |
2395 | | } |
| 2396 | numberKnapsack = 0; |
| 2397 | for (kRow = 0; kRow < numberRows; kRow++) { |
| 2398 | iRow = kRow; |
| 2399 | /* Just obvious one at first |
| 2400 | positive non unit coefficients |
| 2401 | all integer |
| 2402 | positive rowUpper |
| 2403 | for now - linear (but further down in code may use nonlinear) |
| 2404 | column bounds should be tight |
| 2405 | */ |
| 2406 | //double lower = coinModel.getRowLower(iRow); |
| 2407 | double upper = coinModel.getRowUpper(iRow); |
| 2408 | if (upper < 1.0e10) { |
| 2409 | CoinModelLink triple = coinModel.firstInRow(iRow); |
| 2410 | bool possible = true; |
| 2411 | int n = 0; |
| 2412 | int n1 = 0; |
| 2413 | while (triple.column() >= 0) { |
| 2414 | int iColumn = triple.column(); |
| 2415 | const char * el = coinModel.getElementAsString(iRow, iColumn); |
| 2416 | if (!strcmp("Numeric", el)) { |
| 2417 | if (coinModel.columnLower(iColumn) == coinModel.columnUpper(iColumn)) { |
| 2418 | triple = coinModel.next(triple); |
| 2419 | continue; // fixed |
| 2420 | } |
| 2421 | double value = coinModel.getElement(iRow, iColumn); |
| 2422 | if (value < 0.0) { |
| 2423 | possible = false; |
| 2424 | } else { |
| 2425 | n++; |
| 2426 | if (value == 1.0) |
| 2427 | n1++; |
| 2428 | if (coinModel.columnLower(iColumn) < 0.0) |
| 2429 | possible = false; |
| 2430 | if (!coinModel.isInteger(iColumn)) |
| 2431 | possible = false; |
| 2432 | if (whichKnapsack[iColumn] >= 0) |
| 2433 | possible = false; |
| 2434 | } |
| 2435 | } else { |
| 2436 | possible = false; // non linear |
| 2437 | } |
| 2438 | triple = coinModel.next(triple); |
| 2439 | } |
| 2440 | if (n - n1 > 1 && possible) { |
| 2441 | // try |
| 2442 | CoinModelLink triple = coinModel.firstInRow(iRow); |
| 2443 | while (triple.column() >= 0) { |
| 2444 | int iColumn = triple.column(); |
| 2445 | if (coinModel.columnLower(iColumn) != coinModel.columnUpper(iColumn)) |
| 2446 | whichKnapsack[iColumn] = numberKnapsack; |
| 2447 | triple = coinModel.next(triple); |
| 2448 | } |
| 2449 | knapsackRow[numberKnapsack++] = iRow; |
| 2450 | } |
| 2451 | } |
2454 | | } |
2455 | | if (logLevel>0) |
2456 | | printf("%d out of %d candidate rows are possible\n",numberKnapsack,numberRows); |
2457 | | // Check whether we can get rid of nonlinearities |
2458 | | /* mark rows |
2459 | | -2 in knapsack and other variables |
2460 | | -1 not involved |
2461 | | n only in knapsack n |
2462 | | */ |
2463 | | int * markRow = new int [numberRows]; |
2464 | | for (iRow=0;iRow<numberRows;iRow++) |
2465 | | markRow[iRow]=-1; |
2466 | | int canDo=1; // OK and linear |
2467 | | for (iColumn=0;iColumn<numberColumns;iColumn++) { |
2468 | | CoinModelLink triple=coinModel.firstInColumn(iColumn); |
2469 | | int iKnapsack = whichKnapsack[iColumn]; |
2470 | | bool linear=true; |
2471 | | // See if quadratic objective |
2472 | | const char * expr = coinModel.getColumnObjectiveAsString(iColumn); |
2473 | | if (strcmp(expr,"Numeric")) { |
2474 | | linear=false; |
| 2498 | int * markKnapsack = NULL; |
| 2499 | double * coefficient = NULL; |
| 2500 | double * linear = NULL; |
| 2501 | int * whichRow = NULL; |
| 2502 | int * lookupRow = NULL; |
| 2503 | badModel = (canDo == 0); |
| 2504 | if (numberKnapsack && canDo) { |
| 2505 | /* double check - OK if |
| 2506 | no nonlinear |
| 2507 | nonlinear only on columns in knapsack |
| 2508 | nonlinear only on columns in knapsack * ONE other - same for all in knapsack |
| 2509 | AND that is only row connected to knapsack |
| 2510 | (theoretically could split knapsack if two other and small numbers) |
| 2511 | also ONE could be ONE expression - not just a variable |
| 2512 | */ |
| 2513 | int iKnapsack; |
| 2514 | markKnapsack = new int [numberKnapsack]; |
| 2515 | coefficient = new double [numberKnapsack]; |
| 2516 | linear = new double [numberColumns]; |
| 2517 | for (iKnapsack = 0; iKnapsack < numberKnapsack; iKnapsack++) |
| 2518 | markKnapsack[iKnapsack] = -1; |
| 2519 | if (canDo == 2) { |
| 2520 | for (iRow = -1; iRow < numberRows; iRow++) { |
| 2521 | int numberOdd; |
| 2522 | CoinPackedMatrix * row = coinModel.quadraticRow(iRow, linear, numberOdd); |
| 2523 | if (row) { |
| 2524 | // see if valid |
| 2525 | const double * element = row->getElements(); |
| 2526 | const int * column = row->getIndices(); |
| 2527 | const CoinBigIndex * columnStart = row->getVectorStarts(); |
| 2528 | const int * columnLength = row->getVectorLengths(); |
| 2529 | int numberLook = row->getNumCols(); |
| 2530 | for (int i = 0; i < numberLook; i++) { |
| 2531 | int iKnapsack = whichKnapsack[i]; |
| 2532 | if (iKnapsack < 0) { |
| 2533 | // might be able to swap - but for now can't have knapsack in |
| 2534 | for (int j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) { |
| 2535 | int iColumn = column[j]; |
| 2536 | if (whichKnapsack[iColumn] >= 0) { |
| 2537 | canDo = 0; // no good |
| 2538 | badModel = true; |
| 2539 | break; |
| 2540 | } |
| 2541 | } |
| 2542 | } else { |
| 2543 | // OK if in same knapsack - or maybe just one |
| 2544 | int marked = markKnapsack[iKnapsack]; |
| 2545 | for (int j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) { |
| 2546 | int iColumn = column[j]; |
| 2547 | if (whichKnapsack[iColumn] != iKnapsack && whichKnapsack[iColumn] >= 0) { |
| 2548 | canDo = 0; // no good |
| 2549 | badModel = true; |
| 2550 | break; |
| 2551 | } else if (marked == -1) { |
| 2552 | markKnapsack[iKnapsack] = iColumn; |
| 2553 | marked = iColumn; |
| 2554 | coefficient[iKnapsack] = element[j]; |
| 2555 | coinModel.associateElement(coinModel.columnName(iColumn), 1.0); |
| 2556 | } else if (marked != iColumn) { |
| 2557 | badModel = true; |
| 2558 | canDo = 0; // no good |
| 2559 | break; |
| 2560 | } else { |
| 2561 | // could manage with different coefficients - but for now ... |
| 2562 | assert(coefficient[iKnapsack] == element[j]); |
| 2563 | } |
| 2564 | } |
| 2565 | } |
| 2566 | } |
| 2567 | delete row; |
| 2568 | } |
| 2569 | } |
| 2570 | } |
| 2571 | if (canDo) { |
| 2572 | // for any rows which are cuts |
| 2573 | whichRow = new int [numberRows]; |
| 2574 | lookupRow = new int [numberRows]; |
| 2575 | bool someNonlinear = false; |
| 2576 | double maxCoefficient = 1.0; |
| 2577 | for (iKnapsack = 0; iKnapsack < numberKnapsack; iKnapsack++) { |
| 2578 | if (markKnapsack[iKnapsack] >= 0) { |
| 2579 | someNonlinear = true; |
| 2580 | int iColumn = markKnapsack[iKnapsack]; |
| 2581 | maxCoefficient = CoinMax(maxCoefficient, fabs(coefficient[iKnapsack] * coinModel.columnUpper(iColumn))); |
| 2582 | } |
| 2583 | } |
| 2584 | if (someNonlinear) { |
| 2585 | // associate all columns to stop possible error messages |
| 2586 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 2587 | coinModel.associateElement(coinModel.columnName(iColumn), 1.0); |
| 2588 | } |
| 2589 | } |
| 2590 | ClpSimplex tempModel; |
| 2591 | tempModel.loadProblem(coinModel); |
| 2592 | // Create final model - first without knapsacks |
| 2593 | int nCol = 0; |
| 2594 | int nRow = 0; |
| 2595 | for (iRow = 0; iRow < numberRows; iRow++) { |
| 2596 | if (markRow[iRow] < 0) { |
| 2597 | lookupRow[iRow] = nRow; |
| 2598 | whichRow[nRow++] = iRow; |
| 2599 | } else { |
| 2600 | lookupRow[iRow] = -1; |
| 2601 | } |
| 2602 | } |
| 2603 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 2604 | if (whichKnapsack[iColumn] < 0) |
| 2605 | whichColumn[nCol++] = iColumn; |
| 2606 | } |
| 2607 | ClpSimplex finalModelX(&tempModel, nRow, whichRow, nCol, whichColumn, false, false, false); |
| 2608 | OsiClpSolverInterface finalModelY(&finalModelX, true); |
| 2609 | finalModel = finalModelY.clone(); |
| 2610 | finalModelY.releaseClp(); |
| 2611 | // Put back priorities |
| 2612 | const int * priorities = model.priorities(); |
| 2613 | if (priorities) { |
| 2614 | finalModel->findIntegers(false); |
| 2615 | OsiObject ** objects = finalModel->objects(); |
| 2616 | int numberObjects = finalModel->numberObjects(); |
| 2617 | for (int iObj = 0; iObj < numberObjects; iObj++) { |
| 2618 | int iColumn = objects[iObj]->columnNumber(); |
| 2619 | if (iColumn >= 0 && iColumn < nCol) { |
| 2620 | #ifndef NDEBUG |
| 2621 | OsiSimpleInteger * obj = |
| 2622 | dynamic_cast <OsiSimpleInteger *>(objects[iObj]) ; |
| 2623 | #endif |
| 2624 | assert (obj); |
| 2625 | int iPriority = priorities[whichColumn[iColumn]]; |
| 2626 | if (iPriority > 0) |
| 2627 | objects[iObj]->setPriority(iPriority); |
| 2628 | } |
| 2629 | } |
| 2630 | } |
| 2631 | for (iRow = 0; iRow < numberRows; iRow++) { |
| 2632 | whichRow[iRow] = iRow; |
| 2633 | } |
| 2634 | int numberOther = finalModel->getNumCols(); |
| 2635 | int nLargest = 0; |
| 2636 | int nelLargest = 0; |
| 2637 | int nTotal = 0; |
| 2638 | for (iKnapsack = 0; iKnapsack < numberKnapsack; iKnapsack++) { |
| 2639 | iRow = knapsackRow[iKnapsack]; |
| 2640 | int nCreate = maxTotal; |
| 2641 | int nelCreate = coinModel.expandKnapsack(iRow, nCreate, NULL, NULL, NULL, NULL); |
| 2642 | if (nelCreate < 0) |
| 2643 | badModel = true; |
| 2644 | nTotal += nCreate; |
| 2645 | nLargest = CoinMax(nLargest, nCreate); |
| 2646 | nelLargest = CoinMax(nelLargest, nelCreate); |
| 2647 | } |
| 2648 | if (nTotal > maxTotal) |
| 2649 | badModel = true; |
| 2650 | if (!badModel) { |
| 2651 | // Now arrays for building |
| 2652 | nelLargest = CoinMax(nelLargest, nLargest) + 1; |
| 2653 | double * buildObj = new double [nLargest]; |
| 2654 | double * buildElement = new double [nelLargest]; |
| 2655 | int * buildStart = new int[nLargest+1]; |
| 2656 | int * buildRow = new int[nelLargest]; |
| 2657 | // alow for integers in knapsacks |
| 2658 | OsiObject ** object = new OsiObject * [numberKnapsack+nTotal]; |
| 2659 | int nSOS = 0; |
| 2660 | int nObj = numberKnapsack; |
| 2661 | for (iKnapsack = 0; iKnapsack < numberKnapsack; iKnapsack++) { |
| 2662 | knapsackStart[iKnapsack] = finalModel->getNumCols(); |
| 2663 | iRow = knapsackRow[iKnapsack]; |
| 2664 | int nCreate = 10000; |
| 2665 | coinModel.expandKnapsack(iRow, nCreate, buildObj, buildStart, buildRow, buildElement); |
| 2666 | // Redo row numbers |
| 2667 | for (iColumn = 0; iColumn < nCreate; iColumn++) { |
| 2668 | for (int j = buildStart[iColumn]; j < buildStart[iColumn+1]; j++) { |
| 2669 | int jRow = buildRow[j]; |
| 2670 | jRow = lookupRow[jRow]; |
| 2671 | assert (jRow >= 0 && jRow < nRow); |
| 2672 | buildRow[j] = jRow; |
| 2673 | } |
| 2674 | } |
| 2675 | finalModel->addCols(nCreate, buildStart, buildRow, buildElement, NULL, NULL, buildObj); |
| 2676 | int numberFinal = finalModel->getNumCols(); |
| 2677 | for (iColumn = numberOther; iColumn < numberFinal; iColumn++) { |
| 2678 | if (markKnapsack[iKnapsack] < 0) { |
| 2679 | finalModel->setColUpper(iColumn, maxCoefficient); |
| 2680 | finalModel->setInteger(iColumn); |
| 2681 | } else { |
| 2682 | finalModel->setColUpper(iColumn, maxCoefficient + 1.0); |
| 2683 | finalModel->setInteger(iColumn); |
| 2684 | } |
| 2685 | OsiSimpleInteger * sosObject = new OsiSimpleInteger(finalModel, iColumn); |
| 2686 | sosObject->setPriority(1000000); |
| 2687 | object[nObj++] = sosObject; |
| 2688 | buildRow[iColumn-numberOther] = iColumn; |
| 2689 | buildElement[iColumn-numberOther] = 1.0; |
| 2690 | } |
| 2691 | if (markKnapsack[iKnapsack] < 0) { |
| 2692 | // convexity row |
| 2693 | finalModel->addRow(numberFinal - numberOther, buildRow, buildElement, 1.0, 1.0); |
| 2694 | } else { |
| 2695 | int iColumn = markKnapsack[iKnapsack]; |
| 2696 | int n = numberFinal - numberOther; |
| 2697 | buildRow[n] = iColumn; |
| 2698 | buildElement[n++] = -fabs(coefficient[iKnapsack]); |
| 2699 | // convexity row (sort of) |
| 2700 | finalModel->addRow(n, buildRow, buildElement, 0.0, 0.0); |
| 2701 | OsiSOS * sosObject = new OsiSOS(finalModel, n - 1, buildRow, NULL, 1); |
| 2702 | sosObject->setPriority(iKnapsack + SOSPriority); |
| 2703 | // Say not integral even if is (switch off heuristics) |
| 2704 | sosObject->setIntegerValued(false); |
| 2705 | object[nSOS++] = sosObject; |
| 2706 | } |
| 2707 | numberOther = numberFinal; |
| 2708 | } |
| 2709 | finalModel->addObjects(nObj, object); |
| 2710 | for (iKnapsack = 0; iKnapsack < nObj; iKnapsack++) |
| 2711 | delete object[iKnapsack]; |
| 2712 | delete [] object; |
| 2713 | // Can we move any rows to cuts |
| 2714 | const int * cutMarker = coinModel.cutMarker(); |
| 2715 | if (cutMarker && 0) { |
| 2716 | printf("AMPL CUTS OFF until global cuts fixed\n"); |
| 2717 | cutMarker = NULL; |
| 2718 | } |
| 2719 | if (cutMarker) { |
| 2720 | // Row copy |
| 2721 | const CoinPackedMatrix * matrixByRow = finalModel->getMatrixByRow(); |
| 2722 | const double * elementByRow = matrixByRow->getElements(); |
| 2723 | const int * column = matrixByRow->getIndices(); |
| 2724 | const CoinBigIndex * rowStart = matrixByRow->getVectorStarts(); |
| 2725 | const int * rowLength = matrixByRow->getVectorLengths(); |
| 2726 | |
| 2727 | const double * rowLower = finalModel->getRowLower(); |
| 2728 | const double * rowUpper = finalModel->getRowUpper(); |
| 2729 | int nDelete = 0; |
| 2730 | for (iRow = 0; iRow < numberRows; iRow++) { |
| 2731 | if (cutMarker[iRow] && lookupRow[iRow] >= 0) { |
| 2732 | int jRow = lookupRow[iRow]; |
| 2733 | whichRow[nDelete++] = jRow; |
| 2734 | int start = rowStart[jRow]; |
| 2735 | stored.addCut(rowLower[jRow], rowUpper[jRow], |
| 2736 | rowLength[jRow], column + start, elementByRow + start); |
| 2737 | } |
| 2738 | } |
| 2739 | finalModel->deleteRows(nDelete, whichRow); |
| 2740 | } |
| 2741 | knapsackStart[numberKnapsack] = finalModel->getNumCols(); |
| 2742 | delete [] buildObj; |
| 2743 | delete [] buildElement; |
| 2744 | delete [] buildStart; |
| 2745 | delete [] buildRow; |
| 2746 | finalModel->writeMps("full"); |
| 2747 | } |
| 2748 | } |
2491 | | if (!linear) { |
2492 | | if (whichKnapsack[iColumn]<0) { |
2493 | | canDo=0; |
2494 | | break; |
2495 | | } else { |
2496 | | canDo=2; |
2497 | | } |
2498 | | } |
2499 | | } |
2500 | | int * markKnapsack = NULL; |
2501 | | double * coefficient = NULL; |
2502 | | double * linear = NULL; |
2503 | | int * whichRow = NULL; |
2504 | | int * lookupRow = NULL; |
2505 | | badModel=(canDo==0); |
2506 | | if (numberKnapsack&&canDo) { |
2507 | | /* double check - OK if |
2508 | | no nonlinear |
2509 | | nonlinear only on columns in knapsack |
2510 | | nonlinear only on columns in knapsack * ONE other - same for all in knapsack |
2511 | | AND that is only row connected to knapsack |
2512 | | (theoretically could split knapsack if two other and small numbers) |
2513 | | also ONE could be ONE expression - not just a variable |
2514 | | */ |
2515 | | int iKnapsack; |
2516 | | markKnapsack = new int [numberKnapsack]; |
2517 | | coefficient = new double [numberKnapsack]; |
2518 | | linear = new double [numberColumns]; |
2519 | | for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) |
2520 | | markKnapsack[iKnapsack]=-1; |
2521 | | if (canDo==2) { |
2522 | | for (iRow=-1;iRow<numberRows;iRow++) { |
2523 | | int numberOdd; |
2524 | | CoinPackedMatrix * row = coinModel.quadraticRow(iRow,linear,numberOdd); |
2525 | | if (row) { |
2526 | | // see if valid |
2527 | | const double * element = row->getElements(); |
2528 | | const int * column = row->getIndices(); |
2529 | | const CoinBigIndex * columnStart = row->getVectorStarts(); |
2530 | | const int * columnLength = row->getVectorLengths(); |
2531 | | int numberLook = row->getNumCols(); |
2532 | | for (int i=0;i<numberLook;i++) { |
2533 | | int iKnapsack=whichKnapsack[i]; |
2534 | | if (iKnapsack<0) { |
2535 | | // might be able to swap - but for now can't have knapsack in |
2536 | | for (int j=columnStart[i];j<columnStart[i]+columnLength[i];j++) { |
2537 | | int iColumn = column[j]; |
2538 | | if (whichKnapsack[iColumn]>=0) { |
2539 | | canDo=0; // no good |
2540 | | badModel=true; |
2541 | | break; |
2542 | | } |
2543 | | } |
2544 | | } else { |
2545 | | // OK if in same knapsack - or maybe just one |
2546 | | int marked=markKnapsack[iKnapsack]; |
2547 | | for (int j=columnStart[i];j<columnStart[i]+columnLength[i];j++) { |
2548 | | int iColumn = column[j]; |
2549 | | if (whichKnapsack[iColumn]!=iKnapsack&&whichKnapsack[iColumn]>=0) { |
2550 | | canDo=0; // no good |
2551 | | badModel=true; |
2552 | | break; |
2553 | | } else if (marked==-1) { |
2554 | | markKnapsack[iKnapsack]=iColumn; |
2555 | | marked=iColumn; |
2556 | | coefficient[iKnapsack]=element[j]; |
2557 | | coinModel.associateElement(coinModel.columnName(iColumn),1.0); |
2558 | | } else if (marked!=iColumn) { |
2559 | | badModel=true; |
2560 | | canDo=0; // no good |
2561 | | break; |
2562 | | } else { |
2563 | | // could manage with different coefficients - but for now ... |
2564 | | assert(coefficient[iKnapsack]==element[j]); |
2565 | | } |
2566 | | } |
2567 | | } |
2568 | | } |
2569 | | delete row; |
2570 | | } |
2571 | | } |
2572 | | } |
2573 | | if (canDo) { |
2574 | | // for any rows which are cuts |
2575 | | whichRow = new int [numberRows]; |
2576 | | lookupRow = new int [numberRows]; |
2577 | | bool someNonlinear=false; |
2578 | | double maxCoefficient=1.0; |
2579 | | for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) { |
2580 | | if (markKnapsack[iKnapsack]>=0) { |
2581 | | someNonlinear=true; |
2582 | | int iColumn = markKnapsack[iKnapsack]; |
2583 | | maxCoefficient = CoinMax(maxCoefficient,fabs(coefficient[iKnapsack]*coinModel.columnUpper(iColumn))); |
2584 | | } |
2585 | | } |
2586 | | if (someNonlinear) { |
2587 | | // associate all columns to stop possible error messages |
2588 | | for (iColumn=0;iColumn<numberColumns;iColumn++) { |
2589 | | coinModel.associateElement(coinModel.columnName(iColumn),1.0); |
2590 | | } |
2591 | | } |
2592 | | ClpSimplex tempModel; |
2593 | | tempModel.loadProblem(coinModel); |
2594 | | // Create final model - first without knapsacks |
2595 | | int nCol=0; |
2596 | | int nRow=0; |
2597 | | for (iRow=0;iRow<numberRows;iRow++) { |
2598 | | if (markRow[iRow]<0) { |
2599 | | lookupRow[iRow]=nRow; |
2600 | | whichRow[nRow++]=iRow; |
2601 | | } else { |
2602 | | lookupRow[iRow]=-1; |
2603 | | } |
2604 | | } |
2605 | | for (iColumn=0;iColumn<numberColumns;iColumn++) { |
2606 | | if (whichKnapsack[iColumn]<0) |
2607 | | whichColumn[nCol++]=iColumn; |
2608 | | } |
2609 | | ClpSimplex finalModelX(&tempModel,nRow,whichRow,nCol,whichColumn,false,false,false); |
2610 | | OsiClpSolverInterface finalModelY(&finalModelX,true); |
2611 | | finalModel = finalModelY.clone(); |
2612 | | finalModelY.releaseClp(); |
2613 | | // Put back priorities |
2614 | | const int * priorities=model.priorities(); |
2615 | | if (priorities) { |
2616 | | finalModel->findIntegers(false); |
2617 | | OsiObject ** objects = finalModel->objects(); |
2618 | | int numberObjects = finalModel->numberObjects(); |
2619 | | for (int iObj = 0;iObj<numberObjects;iObj++) { |
2620 | | int iColumn = objects[iObj]->columnNumber(); |
2621 | | if (iColumn>=0&&iColumn<nCol) { |
2622 | | #ifndef NDEBUG |
2623 | | OsiSimpleInteger * obj = |
2624 | | dynamic_cast <OsiSimpleInteger *>(objects[iObj]) ; |
2625 | | #endif |
2626 | | assert (obj); |
2627 | | int iPriority = priorities[whichColumn[iColumn]]; |
2628 | | if (iPriority>0) |
2629 | | objects[iObj]->setPriority(iPriority); |
2630 | | } |
2631 | | } |
2632 | | } |
2633 | | for (iRow=0;iRow<numberRows;iRow++) { |
2634 | | whichRow[iRow]=iRow; |
2635 | | } |
2636 | | int numberOther=finalModel->getNumCols(); |
2637 | | int nLargest=0; |
2638 | | int nelLargest=0; |
2639 | | int nTotal=0; |
2640 | | for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) { |
2641 | | iRow = knapsackRow[iKnapsack]; |
2642 | | int nCreate = maxTotal; |
2643 | | int nelCreate=coinModel.expandKnapsack(iRow,nCreate,NULL,NULL,NULL,NULL); |
2644 | | if (nelCreate<0) |
2645 | | badModel=true; |
2646 | | nTotal+=nCreate; |
2647 | | nLargest = CoinMax(nLargest,nCreate); |
2648 | | nelLargest = CoinMax(nelLargest,nelCreate); |
2649 | | } |
2650 | | if (nTotal>maxTotal) |
2651 | | badModel=true; |
2652 | | if (!badModel) { |
2653 | | // Now arrays for building |
2654 | | nelLargest = CoinMax(nelLargest,nLargest)+1; |
2655 | | double * buildObj = new double [nLargest]; |
2656 | | double * buildElement = new double [nelLargest]; |
2657 | | int * buildStart = new int[nLargest+1]; |
2658 | | int * buildRow = new int[nelLargest]; |
2659 | | // alow for integers in knapsacks |
2660 | | OsiObject ** object = new OsiObject * [numberKnapsack+nTotal]; |
2661 | | int nSOS=0; |
2662 | | int nObj=numberKnapsack; |
2663 | | for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) { |
2664 | | knapsackStart[iKnapsack]=finalModel->getNumCols(); |
2665 | | iRow = knapsackRow[iKnapsack]; |
2666 | | int nCreate = 10000; |
2667 | | coinModel.expandKnapsack(iRow,nCreate,buildObj,buildStart,buildRow,buildElement); |
2668 | | // Redo row numbers |
2669 | | for (iColumn=0;iColumn<nCreate;iColumn++) { |
2670 | | for (int j=buildStart[iColumn];j<buildStart[iColumn+1];j++) { |
2671 | | int jRow=buildRow[j]; |
2672 | | jRow=lookupRow[jRow]; |
2673 | | assert (jRow>=0&&jRow<nRow); |
2674 | | buildRow[j]=jRow; |
2675 | | } |
2676 | | } |
2677 | | finalModel->addCols(nCreate,buildStart,buildRow,buildElement,NULL,NULL,buildObj); |
2678 | | int numberFinal=finalModel->getNumCols(); |
2679 | | for (iColumn=numberOther;iColumn<numberFinal;iColumn++) { |
2680 | | if (markKnapsack[iKnapsack]<0) { |
2681 | | finalModel->setColUpper(iColumn,maxCoefficient); |
2682 | | finalModel->setInteger(iColumn); |
2683 | | } else { |
2684 | | finalModel->setColUpper(iColumn,maxCoefficient+1.0); |
2685 | | finalModel->setInteger(iColumn); |
2686 | | } |
2687 | | OsiSimpleInteger * sosObject = new OsiSimpleInteger(finalModel,iColumn); |
2688 | | sosObject->setPriority(1000000); |
2689 | | object[nObj++]=sosObject; |
2690 | | buildRow[iColumn-numberOther]=iColumn; |
2691 | | buildElement[iColumn-numberOther]=1.0; |
2692 | | } |
2693 | | if (markKnapsack[iKnapsack]<0) { |
2694 | | // convexity row |
2695 | | finalModel->addRow(numberFinal-numberOther,buildRow,buildElement,1.0,1.0); |
2696 | | } else { |
2697 | | int iColumn = markKnapsack[iKnapsack]; |
2698 | | int n=numberFinal-numberOther; |
2699 | | buildRow[n]=iColumn; |
2700 | | buildElement[n++]=-fabs(coefficient[iKnapsack]); |
2701 | | // convexity row (sort of) |
2702 | | finalModel->addRow(n,buildRow,buildElement,0.0,0.0); |
2703 | | OsiSOS * sosObject = new OsiSOS(finalModel,n-1,buildRow,NULL,1); |
2704 | | sosObject->setPriority(iKnapsack+SOSPriority); |
2705 | | // Say not integral even if is (switch off heuristics) |
2706 | | sosObject->setIntegerValued(false); |
2707 | | object[nSOS++]=sosObject; |
2708 | | } |
2709 | | numberOther=numberFinal; |
2710 | | } |
2711 | | finalModel->addObjects(nObj,object); |
2712 | | for (iKnapsack=0;iKnapsack<nObj;iKnapsack++) |
2713 | | delete object[iKnapsack]; |
2714 | | delete [] object; |
2715 | | // Can we move any rows to cuts |
2716 | | const int * cutMarker = coinModel.cutMarker(); |
2717 | | if (cutMarker&&0) { |
2718 | | printf("AMPL CUTS OFF until global cuts fixed\n"); |
2719 | | cutMarker=NULL; |
2720 | | } |
2721 | | if (cutMarker) { |
2722 | | // Row copy |
2723 | | const CoinPackedMatrix * matrixByRow = finalModel->getMatrixByRow(); |
2724 | | const double * elementByRow = matrixByRow->getElements(); |
2725 | | const int * column = matrixByRow->getIndices(); |
2726 | | const CoinBigIndex * rowStart = matrixByRow->getVectorStarts(); |
2727 | | const int * rowLength = matrixByRow->getVectorLengths(); |
2728 | | |
2729 | | const double * rowLower = finalModel->getRowLower(); |
2730 | | const double * rowUpper = finalModel->getRowUpper(); |
2731 | | int nDelete=0; |
2732 | | for (iRow=0;iRow<numberRows;iRow++) { |
2733 | | if (cutMarker[iRow]&&lookupRow[iRow]>=0) { |
2734 | | int jRow=lookupRow[iRow]; |
2735 | | whichRow[nDelete++]=jRow; |
2736 | | int start = rowStart[jRow]; |
2737 | | stored.addCut(rowLower[jRow],rowUpper[jRow], |
2738 | | rowLength[jRow],column+start,elementByRow+start); |
2739 | | } |
2740 | | } |
2741 | | finalModel->deleteRows(nDelete,whichRow); |
2742 | | } |
2743 | | knapsackStart[numberKnapsack]=finalModel->getNumCols(); |
2744 | | delete [] buildObj; |
2745 | | delete [] buildElement; |
2746 | | delete [] buildStart; |
2747 | | delete [] buildRow; |
2748 | | finalModel->writeMps("full"); |
2749 | | } |
2750 | | } |
2751 | | } |
2752 | | delete [] whichKnapsack; |
2753 | | delete [] markRow; |
2754 | | delete [] markKnapsack; |
2755 | | delete [] coefficient; |
2756 | | delete [] linear; |
2757 | | delete [] whichRow; |
2758 | | delete [] lookupRow; |
2759 | | delete si; |
2760 | | si=NULL; |
2761 | | if (!badModel&&finalModel) { |
2762 | | finalModel->setDblParam(OsiObjOffset,coinModel.objectiveOffset()); |
2763 | | return finalModel; |
2764 | | } else { |
2765 | | delete finalModel; |
2766 | | printf("can't make knapsacks - did you set fixedPriority (extra1)\n"); |
2767 | | return NULL; |
2768 | | } |