Changeset 1368 for branches/sandbox/Cbc/src/CbcHeuristicRINS.cpp
 Timestamp:
 Dec 5, 2009 4:45:23 PM (10 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

branches/sandbox/Cbc/src/CbcHeuristicRINS.cpp
r1315 r1368 339 339 memset(used_, 0, numberColumns); 340 340 } 341 // Default Constructor 342 CbcHeuristicRENS::CbcHeuristicRENS() 343 : CbcHeuristic() 344 { 345 numberTries_ = 0; 346 whereFrom_ = 256 + 1; 347 } 348 349 // Constructor with model  assumed before cuts 350 351 CbcHeuristicRENS::CbcHeuristicRENS(CbcModel & model) 352 : CbcHeuristic(model) 353 { 354 numberTries_ = 0; 355 whereFrom_ = 256 + 1; 356 } 357 358 // Destructor 359 CbcHeuristicRENS::~CbcHeuristicRENS () 360 { 361 } 362 363 // Clone 364 CbcHeuristic * 365 CbcHeuristicRENS::clone() const 366 { 367 return new CbcHeuristicRENS(*this); 368 } 369 370 // Assignment operator 371 CbcHeuristicRENS & 372 CbcHeuristicRENS::operator=( const CbcHeuristicRENS & rhs) 373 { 374 if (this != &rhs) { 375 CbcHeuristic::operator=(rhs); 376 numberTries_ = rhs.numberTries_; 377 } 378 return *this; 379 } 380 381 // Copy constructor 382 CbcHeuristicRENS::CbcHeuristicRENS(const CbcHeuristicRENS & rhs) 383 : 384 CbcHeuristic(rhs), 385 numberTries_(rhs.numberTries_) 386 { 387 } 388 // Resets stuff if model changes 389 void 390 CbcHeuristicRENS::resetModel(CbcModel * ) 391 { 392 } 393 int 394 CbcHeuristicRENS::solution(double & solutionValue, 395 double * betterSolution) 396 { 397 int returnCode = 0; 398 const double * bestSolution = model_>bestSolution(); 399 if (numberTries_  (when() < 2 && bestSolution)) 400 return 0; 401 numberTries_++; 402 OsiSolverInterface * solver = model_>solver(); 403 404 int numberIntegers = model_>numberIntegers(); 405 const int * integerVariable = model_>integerVariable(); 406 407 const double * currentSolution = solver>getColSolution(); 408 OsiSolverInterface * newSolver = cloneBut(3); // was model_>continuousSolver()>clone(); 409 const double * colLower = newSolver>getColLower(); 410 const double * colUpper = newSolver>getColUpper(); 411 412 double primalTolerance; 413 solver>getDblParam(OsiPrimalTolerance, primalTolerance); 414 415 int i; 416 int numberFixed = 0; 417 int numberTightened = 0; 418 int numberAtBound = 0; 419 int numberColumns = newSolver>getNumCols(); 420 int numberContinuous = numberColumns  numberIntegers; 421 422 for (i = 0; i < numberIntegers; i++) { 423 int iColumn = integerVariable[i]; 424 double value = currentSolution[iColumn]; 425 double lower = colLower[iColumn]; 426 double upper = colUpper[iColumn]; 427 value = CoinMax(value, lower); 428 value = CoinMin(value, upper); 429 #define RENS_FIX_ONLY_LOWER 430 #ifndef RENS_FIX_ONLY_LOWER 431 if (fabs(value  floor(value + 0.5)) < 1.0e8) { 432 value = floor(value + 0.5); 433 if (value == lower  value == upper) 434 numberAtBound++; 435 newSolver>setColLower(iColumn, value); 436 newSolver>setColUpper(iColumn, value); 437 numberFixed++; 438 } else if (colUpper[iColumn]  colLower[iColumn] >= 2.0) { 439 numberTightened++; 440 newSolver>setColLower(iColumn, floor(value)); 441 newSolver>setColUpper(iColumn, ceil(value)); 442 } 443 #else 444 if (fabs(value  floor(value + 0.5)) < 1.0e8 && 445 floor(value + 0.5) == lower) { 446 value = floor(value + 0.5); 447 numberAtBound++; 448 newSolver>setColLower(iColumn, value); 449 newSolver>setColUpper(iColumn, value); 450 numberFixed++; 451 } else if (colUpper[iColumn]  colLower[iColumn] >= 2.0) { 452 numberTightened++; 453 if (fabs(value  floor(value + 0.5)) < 1.0e8) { 454 value = floor(value + 0.5); 455 if (value < upper) { 456 newSolver>setColLower(iColumn, CoinMax(value  1.0, lower)); 457 newSolver>setColUpper(iColumn, CoinMin(value + 1.0, upper)); 458 } else { 459 newSolver>setColLower(iColumn, upper  1.0); 460 } 461 } else { 462 newSolver>setColLower(iColumn, floor(value)); 463 newSolver>setColUpper(iColumn, ceil(value)); 464 } 465 } 466 #endif 467 } 468 if (numberFixed > numberIntegers / 5) { 469 if (numberContinuous > numberIntegers && numberFixed < numberColumns / 5) { 470 #define RENS_FIX_CONTINUOUS 471 #ifdef RENS_FIX_CONTINUOUS 472 const double * colLower = newSolver>getColLower(); 473 //const double * colUpper = newSolver>getColUpper(); 474 int nAtLb = 0; 475 double sumDj = 0.0; 476 const double * dj = newSolver>getReducedCost(); 477 double direction = newSolver>getObjSense(); 478 for (int iColumn = 0; iColumn < numberColumns; iColumn++) { 479 if (!newSolver>isInteger(iColumn)) { 480 double value = currentSolution[iColumn]; 481 if (value < colLower[iColumn] + 1.0e8) { 482 double djValue = dj[iColumn] * direction; 483 nAtLb++; 484 sumDj += djValue; 485 } 486 } 487 } 488 if (nAtLb) { 489 // fix some continuous 490 double * sort = new double[nAtLb]; 491 int * which = new int [nAtLb]; 492 double threshold = CoinMax((0.01 * sumDj) / static_cast<double>(nAtLb), 1.0e6); 493 int nFix2 = 0; 494 for (int iColumn = 0; iColumn < numberColumns; iColumn++) { 495 if (!newSolver>isInteger(iColumn)) { 496 double value = currentSolution[iColumn]; 497 if (value < colLower[iColumn] + 1.0e8) { 498 double djValue = dj[iColumn] * direction; 499 if (djValue > threshold) { 500 sort[nFix2] = djValue; 501 which[nFix2++] = iColumn; 502 } 503 } 504 } 505 } 506 CoinSort_2(sort, sort + nFix2, which); 507 nFix2 = CoinMin(nFix2, (numberColumns  numberFixed) / 2); 508 for (int i = 0; i < nFix2; i++) { 509 int iColumn = which[i]; 510 newSolver>setColUpper(iColumn, colLower[iColumn]); 511 } 512 delete [] sort; 513 delete [] which; 514 #ifdef CLP_INVESTIGATE2 515 printf("%d integers fixed (%d tightened) (%d at bound), and %d continuous fixed at lb\n", 516 numberFixed, numberTightened, numberAtBound, nFix2); 517 #endif 518 } 519 #endif 520 } 521 #ifdef COIN_DEVELOP 522 printf("%d integers fixed and %d tightened\n", numberFixed, numberTightened); 523 #endif 524 returnCode = smallBranchAndBound(newSolver, numberNodes_, betterSolution, solutionValue, 525 model_>getCutoff(), "CbcHeuristicRENS"); 526 if (returnCode < 0) { 527 returnCode = 0; // returned on size 528 #ifdef RENS_FIX_CONTINUOUS 529 if (numberContinuous > numberIntegers && numberFixed >= numberColumns / 5) { 530 const double * colLower = newSolver>getColLower(); 531 //const double * colUpper = newSolver>getColUpper(); 532 int nAtLb = 0; 533 double sumDj = 0.0; 534 const double * dj = newSolver>getReducedCost(); 535 double direction = newSolver>getObjSense(); 536 for (int iColumn = 0; iColumn < numberColumns; iColumn++) { 537 if (!newSolver>isInteger(iColumn)) { 538 double value = currentSolution[iColumn]; 539 if (value < colLower[iColumn] + 1.0e8) { 540 double djValue = dj[iColumn] * direction; 541 nAtLb++; 542 sumDj += djValue; 543 } 544 } 545 } 546 if (nAtLb) { 547 // fix some continuous 548 double * sort = new double[nAtLb]; 549 int * which = new int [nAtLb]; 550 double threshold = CoinMax((0.01 * sumDj) / static_cast<double>(nAtLb), 1.0e6); 551 int nFix2 = 0; 552 for (int iColumn = 0; iColumn < numberColumns; iColumn++) { 553 if (!newSolver>isInteger(iColumn)) { 554 double value = currentSolution[iColumn]; 555 if (value < colLower[iColumn] + 1.0e8) { 556 double djValue = dj[iColumn] * direction; 557 if (djValue > threshold) { 558 sort[nFix2] = djValue; 559 which[nFix2++] = iColumn; 560 } 561 } 562 } 563 } 564 CoinSort_2(sort, sort + nFix2, which); 565 nFix2 = CoinMin(nFix2, (numberColumns  numberFixed) / 2); 566 for (int i = 0; i < nFix2; i++) { 567 int iColumn = which[i]; 568 newSolver>setColUpper(iColumn, colLower[iColumn]); 569 } 570 delete [] sort; 571 delete [] which; 572 #ifdef CLP_INVESTIGATE2 573 printf("%d integers fixed (%d tightened) (%d at bound), and %d continuous fixed at lb\n", 574 numberFixed, numberTightened, numberAtBound, nFix2); 575 #endif 576 } 577 returnCode = smallBranchAndBound(newSolver, numberNodes_, betterSolution, solutionValue, 578 model_>getCutoff(), "CbcHeuristicRENS"); 579 #endif 580 } 581 } 582 //printf("return code %d",returnCode); 583 if ((returnCode&2) != 0) { 584 // could add cut 585 returnCode &= ~2; 586 #ifdef COIN_DEVELOP 587 if (!numberTightened && numberFixed == numberAtBound) 588 printf("could add cut with %d elements\n", numberFixed); 589 #endif 590 } else { 591 //printf("\n"); 592 } 593 } 594 595 delete newSolver; 596 return returnCode; 597 } 598 // update model 599 void CbcHeuristicRENS::setModel(CbcModel * model) 600 { 601 model_ = model; 602 } 603 604 // Default Constructor 605 CbcHeuristicDINS::CbcHeuristicDINS() 606 : CbcHeuristic() 607 { 608 numberSolutions_ = 0; 609 numberSuccesses_ = 0; 610 numberTries_ = 0; 611 howOften_ = 100; 612 decayFactor_ = 0.5; 613 maximumKeepSolutions_ = 5; 614 numberKeptSolutions_ = 0; 615 numberIntegers_ = 1; 616 localSpace_ = 10; 617 values_ = NULL; 618 } 619 620 // Constructor with model  assumed before cuts 621 622 CbcHeuristicDINS::CbcHeuristicDINS(CbcModel & model) 623 : CbcHeuristic(model) 624 { 625 numberSolutions_ = 0; 626 numberSuccesses_ = 0; 627 numberTries_ = 0; 628 howOften_ = 100; 629 decayFactor_ = 0.5; 630 assert(model.solver()); 631 maximumKeepSolutions_ = 5; 632 numberKeptSolutions_ = 0; 633 numberIntegers_ = 1; 634 localSpace_ = 10; 635 values_ = NULL; 636 } 637 638 // Destructor 639 CbcHeuristicDINS::~CbcHeuristicDINS () 640 { 641 for (int i = 0; i < numberKeptSolutions_; i++) 642 delete [] values_[i]; 643 delete [] values_; 644 } 645 646 // Clone 647 CbcHeuristic * 648 CbcHeuristicDINS::clone() const 649 { 650 return new CbcHeuristicDINS(*this); 651 } 652 653 // Assignment operator 654 CbcHeuristicDINS & 655 CbcHeuristicDINS::operator=( const CbcHeuristicDINS & rhs) 656 { 657 if (this != &rhs) { 658 CbcHeuristic::operator=(rhs); 659 numberSolutions_ = rhs.numberSolutions_; 660 howOften_ = rhs.howOften_; 661 numberSuccesses_ = rhs.numberSuccesses_; 662 numberTries_ = rhs.numberTries_; 663 for (int i = 0; i < numberKeptSolutions_; i++) 664 delete [] values_[i]; 665 delete [] values_; 666 maximumKeepSolutions_ = rhs.maximumKeepSolutions_; 667 numberKeptSolutions_ = rhs.numberKeptSolutions_; 668 numberIntegers_ = rhs.numberIntegers_; 669 localSpace_ = rhs.localSpace_; 670 if (model_ && rhs.values_) { 671 assert (numberIntegers_ >= 0); 672 values_ = new int * [maximumKeepSolutions_]; 673 for (int i = 0; i < maximumKeepSolutions_; i++) 674 values_[i] = CoinCopyOfArray(rhs.values_[i], numberIntegers_); 675 } else { 676 values_ = NULL; 677 } 678 } 679 return *this; 680 } 681 682 // Create C++ lines to get to current state 683 void 684 CbcHeuristicDINS::generateCpp( FILE * fp) 685 { 686 CbcHeuristicDINS other; 687 fprintf(fp, "0#include \"CbcHeuristicDINS.hpp\"\n"); 688 fprintf(fp, "3 CbcHeuristicDINS heuristicDINS(*cbcModel);\n"); 689 CbcHeuristic::generateCpp(fp, "heuristicDINS"); 690 if (howOften_ != other.howOften_) 691 fprintf(fp, "3 heuristicDINS.setHowOften(%d);\n", howOften_); 692 else 693 fprintf(fp, "4 heuristicDINS.setHowOften(%d);\n", howOften_); 694 if (maximumKeepSolutions_ != other.maximumKeepSolutions_) 695 fprintf(fp, "3 heuristicDINS.setMaximumKeep(%d);\n", maximumKeepSolutions_); 696 else 697 fprintf(fp, "4 heuristicDINS.setMaximumKeep(%d);\n", maximumKeepSolutions_); 698 fprintf(fp, "3 cbcModel>addHeuristic(&heuristicDINS);\n"); 699 } 700 701 // Copy constructor 702 CbcHeuristicDINS::CbcHeuristicDINS(const CbcHeuristicDINS & rhs) 703 : 704 CbcHeuristic(rhs), 705 numberSolutions_(rhs.numberSolutions_), 706 howOften_(rhs.howOften_), 707 numberSuccesses_(rhs.numberSuccesses_), 708 numberTries_(rhs.numberTries_), 709 maximumKeepSolutions_(rhs.maximumKeepSolutions_), 710 numberKeptSolutions_(rhs.numberKeptSolutions_), 711 numberIntegers_(rhs.numberIntegers_), 712 localSpace_(rhs.localSpace_) 713 { 714 if (model_ && rhs.values_) { 715 assert (numberIntegers_ >= 0); 716 values_ = new int * [maximumKeepSolutions_]; 717 for (int i = 0; i < maximumKeepSolutions_; i++) 718 values_[i] = CoinCopyOfArray(rhs.values_[i], numberIntegers_); 719 } else { 720 values_ = NULL; 721 } 722 } 723 // Resets stuff if model changes 724 void 725 CbcHeuristicDINS::resetModel(CbcModel * ) 726 { 727 //CbcHeuristic::resetModel(model); 728 for (int i = 0; i < numberKeptSolutions_; i++) 729 delete [] values_[i]; 730 delete [] values_; 731 numberKeptSolutions_ = 0; 732 numberIntegers_ = 1; 733 numberSolutions_ = 0; 734 values_ = NULL; 735 } 736 /* 737 First tries setting a variable to better value. If feasible then 738 tries setting others. If not feasible then tries swaps 739 Returns 1 if solution, 0 if not */ 740 int 741 CbcHeuristicDINS::solution(double & solutionValue, 742 double * betterSolution) 743 { 744 numCouldRun_++; 745 int returnCode = 0; 746 const double * bestSolution = model_>bestSolution(); 747 if (!bestSolution) 748 return 0; // No solution found yet 749 if (numberSolutions_ < model_>getSolutionCount()) { 750 // new solution  add info 751 numberSolutions_ = model_>getSolutionCount(); 752 753 int numberIntegers = model_>numberIntegers(); 754 const int * integerVariable = model_>integerVariable(); 755 if (numberIntegers_ < 0) { 756 numberIntegers_ = numberIntegers; 757 assert (!values_); 758 values_ = new int * [maximumKeepSolutions_]; 759 for (int i = 0; i < maximumKeepSolutions_; i++) 760 values_[i] = NULL; 761 } else { 762 assert (numberIntegers == numberIntegers_); 763 } 764 // move solutions (0 will be most recent) 765 { 766 int * temp = values_[maximumKeepSolutions_1]; 767 for (int i = maximumKeepSolutions_  1; i > 0; i) 768 values_[i] = values_[i1]; 769 if (!temp) 770 temp = new int [numberIntegers_]; 771 values_[0] = temp; 772 } 773 int i; 774 for (i = 0; i < numberIntegers; i++) { 775 int iColumn = integerVariable[i]; 776 double value = bestSolution[iColumn]; 777 double nearest = floor(value + 0.5); 778 values_[0][i] = static_cast<int> (nearest); 779 } 780 numberKeptSolutions_ = CoinMin(numberKeptSolutions_ + 1, maximumKeepSolutions_); 781 } 782 int finalReturnCode = 0; 783 if (((model_>getNodeCount() % howOften_) == howOften_ / 2  !model_>getNodeCount()) && (model_>getCurrentPassNumber() == 1  model_>getCurrentPassNumber() == 999999)) { 784 OsiSolverInterface * solver = model_>solver(); 785 786 int numberIntegers = model_>numberIntegers(); 787 const int * integerVariable = model_>integerVariable(); 788 789 const double * currentSolution = solver>getColSolution(); 790 int localSpace = localSpace_; 791 // 0 means finished but no solution, 1 solution, 2 node limit 792 int status = 1; 793 double cutoff = model_>getCutoff(); 794 while (status) { 795 status = 0; 796 OsiSolverInterface * newSolver = cloneBut(3); // was model_>continuousSolver()>clone(); 797 const double * colLower = solver>getColLower(); 798 const double * colUpper = solver>getColUpper(); 799 800 double primalTolerance; 801 solver>getDblParam(OsiPrimalTolerance, primalTolerance); 802 const double * continuousSolution = newSolver>getColSolution(); 803 // Space for added constraint 804 double * element = new double [numberIntegers]; 805 int * column = new int [numberIntegers]; 806 int i; 807 int nFix = 0; 808 int nCouldFix = 0; 809 int nCouldFix2 = 0; 810 int nBound = 0; 811 int nEl = 0; 812 double bias = localSpace; 813 int okSame = numberKeptSolutions_  1; 814 for (i = 0; i < numberIntegers; i++) { 815 int iColumn = integerVariable[i]; 816 const OsiObject * object = model_>object(i); 817 // get original bounds 818 double originalLower; 819 double originalUpper; 820 getIntegerInformation( object, originalLower, originalUpper); 821 double valueInt = bestSolution[iColumn]; 822 if (valueInt < originalLower) { 823 valueInt = originalLower; 824 } else if (valueInt > originalUpper) { 825 valueInt = originalUpper; 826 } 827 int intValue = static_cast<int> (floor(valueInt + 0.5)); 828 double currentValue = currentSolution[iColumn]; 829 double currentLower = colLower[iColumn]; 830 double currentUpper = colUpper[iColumn]; 831 if (fabs(valueInt  currentValue) >= 0.5) { 832 // Rebound 833 nBound++; 834 if (intValue >= currentValue) { 835 currentLower = CoinMax(currentLower, ceil(2 * currentValue  intValue)); 836 currentUpper = intValue; 837 } else { 838 currentLower = intValue; 839 currentUpper = CoinMin(currentUpper, floor(2 * currentValue  intValue)); 840 } 841 newSolver>setColLower(iColumn, currentLower); 842 newSolver>setColUpper(iColumn, currentUpper); 843 } else { 844 // See if can fix 845 bool canFix = false; 846 double continuousValue = continuousSolution[iColumn]; 847 if (fabs(currentValue  valueInt) < 10.0*primalTolerance) { 848 if (currentUpper  currentLower > 1.0) { 849 // General integer variable 850 canFix = true; 851 } else if (fabs(continuousValue  valueInt) < 10.0*primalTolerance) { 852 int nSame = 1; 853 //assert (intValue==values_[0][i]); 854 for (int k = 1; k < numberKeptSolutions_; k++) { 855 if (intValue == values_[k][i]) 856 nSame++; 857 } 858 if (nSame >= okSame) { 859 // can fix 860 canFix = true; 861 } else { 862 nCouldFix++; 863 } 864 } else { 865 nCouldFix2++; 866 } 867 } 868 if (canFix) { 869 newSolver>setColLower(iColumn, intValue); 870 newSolver>setColUpper(iColumn, intValue); 871 nFix++; 872 } else { 873 if (currentUpper  currentLower > 1.0) { 874 // General integer variable 875 currentLower = floor(currentValue); 876 if (intValue >= currentLower && intValue <= currentLower + 1) { 877 newSolver>setColLower(iColumn, currentLower); 878 newSolver>setColUpper(iColumn, currentLower + 1.0); 879 } else { 880 // fix 881 double value; 882 if (intValue < currentLower) 883 value = currentLower; 884 else 885 value = currentLower + 1; 886 newSolver>setColLower(iColumn, value); 887 newSolver>setColUpper(iColumn, value); 888 nFix++; 889 } 890 } else { 891 // 01 (ish) 892 column[nEl] = iColumn; 893 if (intValue == currentLower) { 894 bias += currentLower; 895 element[nEl++] = 1.0; 896 } else if (intValue == currentUpper) { 897 bias += currentUpper; 898 element[nEl++] = 1.0; 899 } else { 900 printf("bad DINS logic\n"); 901 abort(); 902 } 903 } 904 } 905 } 906 } 907 char generalPrint[200]; 908 sprintf(generalPrint, 909 "%d fixed, %d same as cont/int, %d same as int  %d bounded %d in cut\n", 910 nFix, nCouldFix, nCouldFix2, nBound, nEl); 911 model_>messageHandler()>message(CBC_FPUMP2, model_>messages()) 912 << generalPrint 913 << CoinMessageEol; 914 if (nFix > numberIntegers / 10) { 915 #if 0 916 newSolver>initialSolve(); 917 printf("obj %g\n", newSolver>getObjValue()); 918 for (i = 0; i < numberIntegers; i++) { 919 int iColumn = integerVariable[i]; 920 printf("%d new bounds %g %g  solutions %g %g\n", 921 iColumn, newSolver>getColLower()[iColumn], 922 newSolver>getColUpper()[iColumn], 923 bestSolution[iColumn], 924 currentSolution[iColumn]); 925 } 926 #endif 927 if (nEl > 0) 928 newSolver>addRow(nEl, column, element, COIN_DBL_MAX, bias); 929 //printf("%d integers have same value\n",nFix); 930 returnCode = smallBranchAndBound(newSolver, numberNodes_, betterSolution, solutionValue, 931 cutoff, "CbcHeuristicDINS"); 932 if (returnCode < 0) { 933 returnCode = 0; // returned on size 934 status = 0; 935 } else { 936 numRuns_++; 937 if ((returnCode&2) != 0) { 938 // could add cut as complete search 939 returnCode &= ~2; 940 if ((returnCode&1) != 0) { 941 numberSuccesses_++; 942 status = 1; 943 } else { 944 // no solution 945 status = 0; 946 } 947 } else { 948 if ((returnCode&1) != 0) { 949 numberSuccesses_++; 950 status = 1; 951 } else { 952 // no solution but node limit 953 status = 2; 954 if (nEl) 955 localSpace = 5; 956 else 957 localSpace = 1; 958 if (localSpace < 0) 959 status = 0; 960 } 961 } 962 if ((returnCode&1) != 0) { 963 cutoff = CoinMin(cutoff, solutionValue  model_>getCutoffIncrement()); 964 finalReturnCode = 1; 965 } 966 } 967 } 968 delete [] element; 969 delete [] column; 970 delete newSolver; 971 } 972 numberTries_++; 973 if ((numberTries_ % 10) == 0 && numberSuccesses_*3 < numberTries_) 974 howOften_ += static_cast<int> (howOften_ * decayFactor_); 975 } 976 return finalReturnCode; 977 } 978 // update model 979 void CbcHeuristicDINS::setModel(CbcModel * model) 980 { 981 model_ = model; 982 // Get a copy of original matrix 983 assert(model_>solver()); 984 for (int i = 0; i < numberKeptSolutions_; i++) 985 delete [] values_[i]; 986 delete [] values_; 987 numberKeptSolutions_ = 0; 988 numberIntegers_ = 1; 989 numberSolutions_ = 0; 990 values_ = NULL; 991 } 992 993 // Default Constructor 994 CbcHeuristicVND::CbcHeuristicVND() 995 : CbcHeuristic() 996 { 997 numberSolutions_ = 0; 998 numberSuccesses_ = 0; 999 numberTries_ = 0; 1000 lastNode_ = 999999; 1001 howOften_ = 100; 1002 decayFactor_ = 0.5; 1003 baseSolution_ = NULL; 1004 whereFrom_ = 1 + 8 + 255 * 256; 1005 stepSize_ = 0; 1006 k_ = 0; 1007 kmax_ = 0; 1008 nDifferent_ = 0; 1009 } 1010 1011 // Constructor with model  assumed before cuts 1012 1013 CbcHeuristicVND::CbcHeuristicVND(CbcModel & model) 1014 : CbcHeuristic(model) 1015 { 1016 numberSolutions_ = 0; 1017 numberSuccesses_ = 0; 1018 numberTries_ = 0; 1019 lastNode_ = 999999; 1020 howOften_ = 100; 1021 decayFactor_ = 0.5; 1022 assert(model.solver()); 1023 int numberColumns = model.solver()>getNumCols(); 1024 baseSolution_ = new double [numberColumns]; 1025 memset(baseSolution_, 0, numberColumns*sizeof(double)); 1026 whereFrom_ = 1 + 8 + 255 * 256; 1027 stepSize_ = 0; 1028 k_ = 0; 1029 kmax_ = 0; 1030 nDifferent_ = 0; 1031 } 1032 1033 // Destructor 1034 CbcHeuristicVND::~CbcHeuristicVND () 1035 { 1036 delete [] baseSolution_; 1037 } 1038 1039 // Clone 1040 CbcHeuristic * 1041 CbcHeuristicVND::clone() const 1042 { 1043 return new CbcHeuristicVND(*this); 1044 } 1045 1046 // Assignment operator 1047 CbcHeuristicVND & 1048 CbcHeuristicVND::operator=( const CbcHeuristicVND & rhs) 1049 { 1050 if (this != &rhs) { 1051 CbcHeuristic::operator=(rhs); 1052 numberSolutions_ = rhs.numberSolutions_; 1053 howOften_ = rhs.howOften_; 1054 numberSuccesses_ = rhs.numberSuccesses_; 1055 numberTries_ = rhs.numberTries_; 1056 lastNode_ = rhs.lastNode_; 1057 delete [] baseSolution_; 1058 if (model_ && rhs.baseSolution_) { 1059 int numberColumns = model_>solver()>getNumCols(); 1060 baseSolution_ = new double [numberColumns]; 1061 memcpy(baseSolution_, rhs.baseSolution_, numberColumns*sizeof(double)); 1062 } else { 1063 baseSolution_ = NULL; 1064 } 1065 stepSize_ = rhs.stepSize_; 1066 k_ = rhs.k_; 1067 kmax_ = rhs.kmax_; 1068 nDifferent_ = rhs.nDifferent_; 1069 } 1070 return *this; 1071 } 1072 1073 // Create C++ lines to get to current state 1074 void 1075 CbcHeuristicVND::generateCpp( FILE * fp) 1076 { 1077 CbcHeuristicVND other; 1078 fprintf(fp, "0#include \"CbcHeuristicVND.hpp\"\n"); 1079 fprintf(fp, "3 CbcHeuristicVND heuristicVND(*cbcModel);\n"); 1080 CbcHeuristic::generateCpp(fp, "heuristicVND"); 1081 if (howOften_ != other.howOften_) 1082 fprintf(fp, "3 heuristicVND.setHowOften(%d);\n", howOften_); 1083 else 1084 fprintf(fp, "4 heuristicVND.setHowOften(%d);\n", howOften_); 1085 fprintf(fp, "3 cbcModel>addHeuristic(&heuristicVND);\n"); 1086 } 1087 1088 // Copy constructor 1089 CbcHeuristicVND::CbcHeuristicVND(const CbcHeuristicVND & rhs) 1090 : 1091 CbcHeuristic(rhs), 1092 numberSolutions_(rhs.numberSolutions_), 1093 howOften_(rhs.howOften_), 1094 numberSuccesses_(rhs.numberSuccesses_), 1095 numberTries_(rhs.numberTries_), 1096 lastNode_(rhs.lastNode_) 1097 { 1098 if (model_ && rhs.baseSolution_) { 1099 int numberColumns = model_>solver()>getNumCols(); 1100 baseSolution_ = new double [numberColumns]; 1101 memcpy(baseSolution_, rhs.baseSolution_, numberColumns*sizeof(double)); 1102 } else { 1103 baseSolution_ = NULL; 1104 } 1105 stepSize_ = rhs.stepSize_; 1106 k_ = rhs.k_; 1107 kmax_ = rhs.kmax_; 1108 nDifferent_ = rhs.nDifferent_; 1109 } 1110 // Resets stuff if model changes 1111 void 1112 CbcHeuristicVND::resetModel(CbcModel * /*model*/) 1113 { 1114 //CbcHeuristic::resetModel(model); 1115 delete [] baseSolution_; 1116 if (model_ && baseSolution_) { 1117 int numberColumns = model_>solver()>getNumCols(); 1118 baseSolution_ = new double [numberColumns]; 1119 memset(baseSolution_, 0, numberColumns*sizeof(double)); 1120 } else { 1121 baseSolution_ = NULL; 1122 } 1123 } 1124 /* 1125 First tries setting a variable to better value. If feasible then 1126 tries setting others. If not feasible then tries swaps 1127 Returns 1 if solution, 0 if not */ 1128 int 1129 CbcHeuristicVND::solution(double & solutionValue, 1130 double * betterSolution) 1131 { 1132 numCouldRun_++; 1133 int returnCode = 0; 1134 const double * bestSolution = model_>bestSolution(); 1135 if (!bestSolution) 1136 return 0; // No solution found yet 1137 if (numberSolutions_ < model_>getSolutionCount()) { 1138 // new solution  add info 1139 numberSolutions_ = model_>getSolutionCount(); 1140 1141 int numberIntegers = model_>numberIntegers(); 1142 const int * integerVariable = model_>integerVariable(); 1143 1144 int i; 1145 for (i = 0; i < numberIntegers; i++) { 1146 int iColumn = integerVariable[i]; 1147 const OsiObject * object = model_>object(i); 1148 // get original bounds 1149 double originalLower; 1150 double originalUpper; 1151 getIntegerInformation( object, originalLower, originalUpper); 1152 double value = bestSolution[iColumn]; 1153 if (value < originalLower) { 1154 value = originalLower; 1155 } else if (value > originalUpper) { 1156 value = originalUpper; 1157 } 1158 } 1159 } 1160 int numberNodes = model_>getNodeCount(); 1161 if (howOften_ == 100) { 1162 if (numberNodes < lastNode_ + 12) 1163 return 0; 1164 // Do at 50 and 100 1165 if ((numberNodes > 40 && numberNodes <= 50)  (numberNodes > 90 && numberNodes < 100)) 1166 numberNodes = howOften_; 1167 } 1168 if ((numberNodes % howOften_) == 0 && (model_>getCurrentPassNumber() == 1  1169 model_>getCurrentPassNumber() == 999999)) { 1170 lastNode_ = model_>getNodeCount(); 1171 OsiSolverInterface * solver = model_>solver(); 1172 1173 int numberIntegers = model_>numberIntegers(); 1174 const int * integerVariable = model_>integerVariable(); 1175 1176 const double * currentSolution = solver>getColSolution(); 1177 OsiSolverInterface * newSolver = cloneBut(3); // was model_>continuousSolver()>clone(); 1178 //const double * colLower = newSolver>getColLower(); 1179 //const double * colUpper = newSolver>getColUpper(); 1180 1181 double primalTolerance; 1182 solver>getDblParam(OsiPrimalTolerance, primalTolerance); 1183 1184 // Sort on distance 1185 double * distance = new double [numberIntegers]; 1186 int * which = new int [numberIntegers]; 1187 1188 int i; 1189 int nFix = 0; 1190 double tolerance = 10.0 * primalTolerance; 1191 for (i = 0; i < numberIntegers; i++) { 1192 int iColumn = integerVariable[i]; 1193 const OsiObject * object = model_>object(i); 1194 // get original bounds 1195 double originalLower; 1196 double originalUpper; 1197 getIntegerInformation( object, originalLower, originalUpper); 1198 double valueInt = bestSolution[iColumn]; 1199 if (valueInt < originalLower) { 1200 valueInt = originalLower; 1201 } else if (valueInt > originalUpper) { 1202 valueInt = originalUpper; 1203 } 1204 baseSolution_[iColumn] = currentSolution[iColumn]; 1205 distance[i] = fabs(currentSolution[iColumn]  valueInt); 1206 which[i] = i; 1207 if (fabs(currentSolution[iColumn]  valueInt) < tolerance) 1208 nFix++; 1209 } 1210 CoinSort_2(distance, distance + numberIntegers, which); 1211 nDifferent_ = numberIntegers  nFix; 1212 stepSize_ = nDifferent_ / 10; 1213 k_ = stepSize_; 1214 //nFix = numberIntegersstepSize_; 1215 for (i = 0; i < nFix; i++) { 1216 int j = which[i]; 1217 int iColumn = integerVariable[j]; 1218 const OsiObject * object = model_>object(i); 1219 // get original bounds 1220 double originalLower; 1221 double originalUpper; 1222 getIntegerInformation( object, originalLower, originalUpper); 1223 double valueInt = bestSolution[iColumn]; 1224 if (valueInt < originalLower) { 1225 valueInt = originalLower; 1226 } else if (valueInt > originalUpper) { 1227 valueInt = originalUpper; 1228 } 1229 double nearest = floor(valueInt + 0.5); 1230 newSolver>setColLower(iColumn, nearest); 1231 newSolver>setColUpper(iColumn, nearest); 1232 } 1233 delete [] distance; 1234 delete [] which; 1235 if (nFix > numberIntegers / 5) { 1236 //printf("%d integers have samish value\n",nFix); 1237 returnCode = smallBranchAndBound(newSolver, numberNodes_, betterSolution, solutionValue, 1238 model_>getCutoff(), "CbcHeuristicVND"); 1239 if (returnCode < 0) 1240 returnCode = 0; // returned on size 1241 else 1242 numRuns_++; 1243 if ((returnCode&1) != 0) 1244 numberSuccesses_++; 1245 //printf("return code %d",returnCode); 1246 if ((returnCode&2) != 0) { 1247 // could add cut 1248 returnCode &= ~2; 1249 //printf("could add cut with %d elements (if all 01)\n",nFix); 1250 } else { 1251 //printf("\n"); 1252 } 1253 numberTries_++; 1254 if ((numberTries_ % 10) == 0 && numberSuccesses_*3 < numberTries_) 1255 howOften_ += static_cast<int> (howOften_ * decayFactor_); 1256 } 1257 1258 delete newSolver; 1259 } 1260 return returnCode; 1261 } 1262 // update model 1263 void CbcHeuristicVND::setModel(CbcModel * model) 1264 { 1265 model_ = model; 1266 // Get a copy of original matrix 1267 assert(model_>solver()); 1268 delete [] baseSolution_; 1269 int numberColumns = model>solver()>getNumCols(); 1270 baseSolution_ = new double [numberColumns]; 1271 memset(baseSolution_, 0, numberColumns*sizeof(double)); 1272 } 1273 1274 341 342 343
Note: See TracChangeset
for help on using the changeset viewer.