58 | | columnArray_[1]->clear(); |
59 | | #endif |
60 | | // long enough for rows+columns |
61 | | int * backPivot = new int [numberRows_+numberColumns_]; |
62 | | int i; |
63 | | for ( i = 0; i < numberRows_ + numberColumns_; i++) { |
64 | | backPivot[i] = -1; |
65 | | } |
66 | | for (i = 0; i < numberRows_; i++) { |
67 | | int iSequence = pivotVariable_[i]; |
68 | | backPivot[iSequence] = i; |
69 | | } |
70 | | // dualTolerance may be zero if from CBC. In fact use that fact |
71 | | bool inCBC = !dualTolerance_; |
72 | | if (inCBC) |
73 | | assert (integerType_); |
74 | | dualTolerance_ = dblParam_[ClpDualTolerance]; |
75 | | double * arrayX = rowArray_[0]->denseVector(); |
76 | | for ( i = 0; i < numberCheck; i++) { |
77 | | rowArray_[0]->clear(); |
78 | | //rowArray_[0]->checkClear(); |
79 | | //rowArray_[1]->checkClear(); |
80 | | //columnArray_[1]->checkClear(); |
81 | | columnArray_[0]->clear(); |
82 | | //columnArray_[0]->checkClear(); |
83 | | int iSequence = which[i]; |
84 | | if (iSequence < 0) { |
85 | | costIncreased[i] = 0.0; |
86 | | sequenceIncreased[i] = -1; |
87 | | costDecreased[i] = 0.0; |
88 | | sequenceDecreased[i] = -1; |
89 | | continue; |
90 | | } |
91 | | double costIncrease = COIN_DBL_MAX; |
92 | | double costDecrease = COIN_DBL_MAX; |
93 | | int sequenceIncrease = -1; |
94 | | int sequenceDecrease = -1; |
95 | | if (valueIncrease) { |
96 | | assert (valueDecrease); |
97 | | valueIncrease[i] = iSequence < numberColumns_ ? columnActivity_[iSequence] : rowActivity_[iSequence-numberColumns_]; |
98 | | valueDecrease[i] = valueIncrease[i]; |
99 | | } |
| 58 | columnArray_[1]->clear(); |
| 59 | #endif |
| 60 | // long enough for rows+columns |
| 61 | int *backPivot = new int[numberRows_ + numberColumns_]; |
| 62 | int i; |
| 63 | for (i = 0; i < numberRows_ + numberColumns_; i++) { |
| 64 | backPivot[i] = -1; |
| 65 | } |
| 66 | for (i = 0; i < numberRows_; i++) { |
| 67 | int iSequence = pivotVariable_[i]; |
| 68 | backPivot[iSequence] = i; |
| 69 | } |
| 70 | // dualTolerance may be zero if from CBC. In fact use that fact |
| 71 | bool inCBC = !dualTolerance_; |
| 72 | if (inCBC) |
| 73 | assert(integerType_); |
| 74 | dualTolerance_ = dblParam_[ClpDualTolerance]; |
| 75 | double *arrayX = rowArray_[0]->denseVector(); |
| 76 | for (i = 0; i < numberCheck; i++) { |
| 77 | rowArray_[0]->clear(); |
| 78 | //rowArray_[0]->checkClear(); |
| 79 | //rowArray_[1]->checkClear(); |
| 80 | //columnArray_[1]->checkClear(); |
| 81 | columnArray_[0]->clear(); |
| 82 | //columnArray_[0]->checkClear(); |
| 83 | int iSequence = which[i]; |
| 84 | if (iSequence < 0) { |
| 85 | costIncreased[i] = 0.0; |
| 86 | sequenceIncreased[i] = -1; |
| 87 | costDecreased[i] = 0.0; |
| 88 | sequenceDecreased[i] = -1; |
| 89 | continue; |
| 90 | } |
| 91 | double costIncrease = COIN_DBL_MAX; |
| 92 | double costDecrease = COIN_DBL_MAX; |
| 93 | int sequenceIncrease = -1; |
| 94 | int sequenceDecrease = -1; |
| 95 | if (valueIncrease) { |
| 96 | assert(valueDecrease); |
| 97 | valueIncrease[i] = iSequence < numberColumns_ ? columnActivity_[iSequence] : rowActivity_[iSequence - numberColumns_]; |
| 98 | valueDecrease[i] = valueIncrease[i]; |
| 99 | } |
250 | | abort(); |
251 | | } |
252 | | } |
253 | | rowArray_[0]->clear(); |
254 | | //rowArray_[1]->clear(); |
255 | | //columnArray_[1]->clear(); |
256 | | columnArray_[0]->clear(); |
257 | | delete [] backPivot; |
258 | | if (!optimizationDirection_) |
259 | | printf("*** ????? Ranging with zero optimization costs\n"); |
| 183 | costDecrease = 0.0; |
| 184 | } |
| 185 | } |
| 186 | costIncrease *= scale2; |
| 187 | costDecrease *= scale2; |
| 188 | } |
| 189 | } break; |
| 190 | case isFixed: |
| 191 | break; |
| 192 | case isFree: |
| 193 | case superBasic: |
| 194 | costIncrease = 0.0; |
| 195 | costDecrease = 0.0; |
| 196 | sequenceIncrease = iSequence; |
| 197 | sequenceDecrease = iSequence; |
| 198 | break; |
| 199 | case atUpperBound: |
| 200 | costIncrease = CoinMax(0.0, -dj_[iSequence]); |
| 201 | sequenceIncrease = iSequence; |
| 202 | if (valueIncrease) |
| 203 | valueIncrease[i] = primalRanging1(iSequence, iSequence); |
| 204 | break; |
| 205 | case atLowerBound: |
| 206 | costDecrease = CoinMax(0.0, dj_[iSequence]); |
| 207 | sequenceDecrease = iSequence; |
| 208 | if (valueIncrease) |
| 209 | valueDecrease[i] = primalRanging1(iSequence, iSequence); |
| 210 | break; |
| 211 | } |
| 212 | double scaleFactor; |
| 213 | if (rowScale_) { |
| 214 | if (iSequence < numberColumns_) |
| 215 | scaleFactor = 1.0 / (objectiveScale_ * columnScale_[iSequence]); |
| 216 | else |
| 217 | scaleFactor = rowScale_[iSequence - numberColumns_] / objectiveScale_; |
| 218 | } else { |
| 219 | scaleFactor = 1.0 / objectiveScale_; |
| 220 | } |
| 221 | if (costIncrease < 1.0e30) |
| 222 | costIncrease *= scaleFactor; |
| 223 | if (costDecrease < 1.0e30) |
| 224 | costDecrease *= scaleFactor; |
| 225 | if (optimizationDirection_ == 1.0) { |
| 226 | costIncreased[i] = costIncrease; |
| 227 | sequenceIncreased[i] = sequenceIncrease; |
| 228 | costDecreased[i] = costDecrease; |
| 229 | sequenceDecreased[i] = sequenceDecrease; |
| 230 | } else if (optimizationDirection_ == -1.0) { |
| 231 | costIncreased[i] = costDecrease; |
| 232 | sequenceIncreased[i] = sequenceDecrease; |
| 233 | costDecreased[i] = costIncrease; |
| 234 | sequenceDecreased[i] = sequenceIncrease; |
| 235 | if (valueIncrease) { |
| 236 | double temp = valueIncrease[i]; |
| 237 | valueIncrease[i] = valueDecrease[i]; |
| 238 | valueDecrease[i] = temp; |
| 239 | } |
| 240 | } else if (optimizationDirection_ == 0.0) { |
| 241 | // !!!!!! ??? |
| 242 | costIncreased[i] = COIN_DBL_MAX; |
| 243 | sequenceIncreased[i] = -1; |
| 244 | costDecreased[i] = COIN_DBL_MAX; |
| 245 | sequenceDecreased[i] = -1; |
| 246 | } else { |
| 247 | abort(); |
| 248 | } |
| 249 | } |
| 250 | rowArray_[0]->clear(); |
| 251 | //rowArray_[1]->clear(); |
| 252 | //columnArray_[1]->clear(); |
| 253 | columnArray_[0]->clear(); |
| 254 | delete[] backPivot; |
| 255 | if (!optimizationDirection_) |
| 256 | printf("*** ????? Ranging with zero optimization costs\n"); |
422 | | unpack(rowArray_[1], iSequence); |
423 | | #endif |
424 | | factorization_->updateColumn(rowArray_[2], rowArray_[1]); |
425 | | // Get extra rows |
426 | | matrix_->extendUpdated(this, rowArray_[1], 0); |
427 | | // do ratio test |
428 | | checkPrimalRatios(rowArray_[1], 1); |
429 | | if (pivotRow_ >= 0) { |
430 | | valueIncrease = theta_; |
431 | | sequenceIncrease = pivotVariable_[pivotRow_]; |
432 | | } |
433 | | checkPrimalRatios(rowArray_[1], -1); |
434 | | if (pivotRow_ >= 0) { |
435 | | valueDecrease = theta_; |
436 | | sequenceDecrease = pivotVariable_[pivotRow_]; |
437 | | } |
438 | | rowArray_[1]->clear(); |
439 | | } |
440 | | break; |
441 | | } |
442 | | double scaleFactor; |
443 | | if (rowScale_) { |
444 | | if (iSequence < numberColumns_) |
445 | | scaleFactor = columnScale_[iSequence] / rhsScale_; |
446 | | else |
447 | | scaleFactor = 1.0 / (rowScale_[iSequence-numberColumns_] * rhsScale_); |
448 | | } else { |
449 | | scaleFactor = 1.0 / rhsScale_; |
450 | | } |
451 | | if (valueIncrease < 1.0e30) |
452 | | valueIncrease *= scaleFactor; |
453 | | else |
454 | | valueIncrease = COIN_DBL_MAX; |
455 | | if (valueDecrease < 1.0e30) |
456 | | valueDecrease *= scaleFactor; |
457 | | else |
458 | | valueDecrease = COIN_DBL_MAX; |
459 | | valueIncreased[i] = valueIncrease; |
460 | | sequenceIncreased[i] = sequenceIncrease; |
461 | | valueDecreased[i] = valueDecrease; |
462 | | sequenceDecreased[i] = sequenceDecrease; |
463 | | } |
| 417 | unpack(rowArray_[1], iSequence); |
| 418 | #endif |
| 419 | factorization_->updateColumn(rowArray_[2], rowArray_[1]); |
| 420 | // Get extra rows |
| 421 | matrix_->extendUpdated(this, rowArray_[1], 0); |
| 422 | // do ratio test |
| 423 | checkPrimalRatios(rowArray_[1], 1); |
| 424 | if (pivotRow_ >= 0) { |
| 425 | valueIncrease = theta_; |
| 426 | sequenceIncrease = pivotVariable_[pivotRow_]; |
| 427 | } |
| 428 | checkPrimalRatios(rowArray_[1], -1); |
| 429 | if (pivotRow_ >= 0) { |
| 430 | valueDecrease = theta_; |
| 431 | sequenceDecrease = pivotVariable_[pivotRow_]; |
| 432 | } |
| 433 | rowArray_[1]->clear(); |
| 434 | } break; |
| 435 | } |
| 436 | double scaleFactor; |
| 437 | if (rowScale_) { |
| 438 | if (iSequence < numberColumns_) |
| 439 | scaleFactor = columnScale_[iSequence] / rhsScale_; |
| 440 | else |
| 441 | scaleFactor = 1.0 / (rowScale_[iSequence - numberColumns_] * rhsScale_); |
| 442 | } else { |
| 443 | scaleFactor = 1.0 / rhsScale_; |
| 444 | } |
| 445 | if (valueIncrease < 1.0e30) |
| 446 | valueIncrease *= scaleFactor; |
| 447 | else |
| 448 | valueIncrease = COIN_DBL_MAX; |
| 449 | if (valueDecrease < 1.0e30) |
| 450 | valueDecrease *= scaleFactor; |
| 451 | else |
| 452 | valueDecrease = COIN_DBL_MAX; |
| 453 | valueIncreased[i] = valueIncrease; |
| 454 | sequenceIncreased[i] = sequenceIncrease; |
| 455 | valueDecreased[i] = valueDecrease; |
| 456 | sequenceDecreased[i] = sequenceDecrease; |
| 457 | } |
666 | | // Set locale so won't get , instead of . |
667 | | char * saveLocale = strdup(setlocale(LC_ALL,NULL)); |
668 | | setlocale(LC_ALL,"C"); |
669 | | if (strcmp(strParam_[ClpProbName].c_str(), "") == 0) { |
670 | | fprintf(fp, "NAME BLANK "); |
671 | | } else { |
672 | | fprintf(fp, "NAME %s ", strParam_[ClpProbName].c_str()); |
673 | | } |
674 | | if (formatType >= 2) |
675 | | fprintf(fp, "FREEIEEE"); |
676 | | else if (writeValues) |
677 | | fprintf(fp, "VALUES"); |
678 | | // finish off name |
679 | | fprintf(fp, "\n"); |
680 | | int iRow = 0; |
681 | | for(int iColumn = 0; iColumn < numberColumns_; iColumn++) { |
682 | | bool printit = false; |
683 | | if( getColumnStatus(iColumn) == ClpSimplex::basic) { |
684 | | printit = true; |
685 | | // Find non basic row |
686 | | for(; iRow < numberRows_; iRow++) { |
687 | | if (getRowStatus(iRow) != ClpSimplex::basic) |
688 | | break; |
689 | | } |
690 | | if (lengthNames_) { |
691 | | if (iRow != numberRows_) { |
692 | | fprintf(fp, " %s %-8s %s", |
693 | | getRowStatus(iRow) == ClpSimplex::atUpperBound ? "XU" : "XL", |
694 | | columnNames_[iColumn].c_str(), |
695 | | rowNames_[iRow].c_str()); |
696 | | iRow++; |
697 | | } else { |
698 | | // Allow for too many basics! |
699 | | fprintf(fp, " BS %-8s ", |
700 | | columnNames_[iColumn].c_str()); |
701 | | // Dummy row name if values |
702 | | if (writeValues) |
703 | | fprintf(fp, " _dummy_"); |
704 | | } |
705 | | } else { |
706 | | // no names |
707 | | if (iRow != numberRows_) { |
708 | | fprintf(fp, " %s C%7.7d R%7.7d", |
709 | | getRowStatus(iRow) == ClpSimplex::atUpperBound ? "XU" : "XL", |
710 | | iColumn, iRow); |
711 | | iRow++; |
712 | | } else { |
713 | | // Allow for too many basics! |
714 | | fprintf(fp, " BS C%7.7d", iColumn); |
715 | | // Dummy row name if values |
716 | | if (writeValues) |
717 | | fprintf(fp, " _dummy_"); |
718 | | } |
719 | | } |
720 | | } else { |
721 | | if( getColumnStatus(iColumn) == ClpSimplex::atUpperBound) { |
722 | | printit = true; |
723 | | if (lengthNames_) |
724 | | fprintf(fp, " UL %s", columnNames_[iColumn].c_str()); |
725 | | else |
726 | | fprintf(fp, " UL C%7.7d", iColumn); |
727 | | // Dummy row name if values |
728 | | if (writeValues) |
729 | | fprintf(fp, " _dummy_"); |
730 | | } else if( (getColumnStatus(iColumn) == ClpSimplex::superBasic|| |
731 | | getColumnStatus(iColumn) == ClpSimplex::isFree)&& |
732 | | writeValues) { |
733 | | printit = true; |
734 | | if (lengthNames_) |
735 | | fprintf(fp, " BS %s", columnNames_[iColumn].c_str()); |
736 | | else |
737 | | fprintf(fp, " BS C%7.7d", iColumn); |
738 | | // Dummy row name if values |
739 | | if (writeValues) |
740 | | fprintf(fp, " _dummy_"); |
741 | | } |
742 | | } |
743 | | if (printit && writeValues) { |
744 | | // add value |
745 | | CoinConvertDouble(0, formatType, columnActivity_[iColumn], number); |
746 | | fprintf(fp, " %s", number); |
747 | | } |
748 | | if (printit) |
749 | | fprintf(fp, "\n"); |
750 | | } |
751 | | fprintf(fp, "ENDATA\n"); |
752 | | fclose(fp); |
753 | | setlocale(LC_ALL,saveLocale); |
754 | | free(saveLocale); |
755 | | return 0; |
| 658 | // Set locale so won't get , instead of . |
| 659 | char *saveLocale = strdup(setlocale(LC_ALL, NULL)); |
| 660 | setlocale(LC_ALL, "C"); |
| 661 | if (strcmp(strParam_[ClpProbName].c_str(), "") == 0) { |
| 662 | fprintf(fp, "NAME BLANK "); |
| 663 | } else { |
| 664 | fprintf(fp, "NAME %s ", strParam_[ClpProbName].c_str()); |
| 665 | } |
| 666 | if (formatType >= 2) |
| 667 | fprintf(fp, "FREEIEEE"); |
| 668 | else if (writeValues) |
| 669 | fprintf(fp, "VALUES"); |
| 670 | // finish off name |
| 671 | fprintf(fp, "\n"); |
| 672 | int iRow = 0; |
| 673 | for (int iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 674 | bool printit = false; |
| 675 | if (getColumnStatus(iColumn) == ClpSimplex::basic) { |
| 676 | printit = true; |
| 677 | // Find non basic row |
| 678 | for (; iRow < numberRows_; iRow++) { |
| 679 | if (getRowStatus(iRow) != ClpSimplex::basic) |
| 680 | break; |
| 681 | } |
| 682 | if (lengthNames_) { |
| 683 | if (iRow != numberRows_) { |
| 684 | fprintf(fp, " %s %-8s %s", |
| 685 | getRowStatus(iRow) == ClpSimplex::atUpperBound ? "XU" : "XL", |
| 686 | columnNames_[iColumn].c_str(), |
| 687 | rowNames_[iRow].c_str()); |
| 688 | iRow++; |
| 689 | } else { |
| 690 | // Allow for too many basics! |
| 691 | fprintf(fp, " BS %-8s ", |
| 692 | columnNames_[iColumn].c_str()); |
| 693 | // Dummy row name if values |
| 694 | if (writeValues) |
| 695 | fprintf(fp, " _dummy_"); |
| 696 | } |
| 697 | } else { |
| 698 | // no names |
| 699 | if (iRow != numberRows_) { |
| 700 | fprintf(fp, " %s C%7.7d R%7.7d", |
| 701 | getRowStatus(iRow) == ClpSimplex::atUpperBound ? "XU" : "XL", |
| 702 | iColumn, iRow); |
| 703 | iRow++; |
| 704 | } else { |
| 705 | // Allow for too many basics! |
| 706 | fprintf(fp, " BS C%7.7d", iColumn); |
| 707 | // Dummy row name if values |
| 708 | if (writeValues) |
| 709 | fprintf(fp, " _dummy_"); |
| 710 | } |
| 711 | } |
| 712 | } else { |
| 713 | if (getColumnStatus(iColumn) == ClpSimplex::atUpperBound) { |
| 714 | printit = true; |
| 715 | if (lengthNames_) |
| 716 | fprintf(fp, " UL %s", columnNames_[iColumn].c_str()); |
| 717 | else |
| 718 | fprintf(fp, " UL C%7.7d", iColumn); |
| 719 | // Dummy row name if values |
| 720 | if (writeValues) |
| 721 | fprintf(fp, " _dummy_"); |
| 722 | } else if ((getColumnStatus(iColumn) == ClpSimplex::superBasic || getColumnStatus(iColumn) == ClpSimplex::isFree) && writeValues) { |
| 723 | printit = true; |
| 724 | if (lengthNames_) |
| 725 | fprintf(fp, " BS %s", columnNames_[iColumn].c_str()); |
| 726 | else |
| 727 | fprintf(fp, " BS C%7.7d", iColumn); |
| 728 | // Dummy row name if values |
| 729 | if (writeValues) |
| 730 | fprintf(fp, " _dummy_"); |
| 731 | } |
| 732 | } |
| 733 | if (printit && writeValues) { |
| 734 | // add value |
| 735 | CoinConvertDouble(0, formatType, columnActivity_[iColumn], number); |
| 736 | fprintf(fp, " %s", number); |
| 737 | } |
| 738 | if (printit) |
| 739 | fprintf(fp, "\n"); |
| 740 | } |
| 741 | fprintf(fp, "ENDATA\n"); |
| 742 | fclose(fp); |
| 743 | setlocale(LC_ALL, saveLocale); |
| 744 | free(saveLocale); |
| 745 | return 0; |
761 | | int status = 0; |
762 | | if (strcmp(fileName, "-") != 0 && strcmp(fileName, "stdin") != 0) { |
763 | | FILE *fp = fopen(fileName, "r"); |
764 | | if (fp) { |
765 | | // can open - lets go for it |
766 | | fclose(fp); |
767 | | } else { |
768 | | handler_->message(CLP_UNABLE_OPEN, messages_) |
769 | | << fileName << CoinMessageEol; |
770 | | return -1; |
771 | | } |
772 | | } |
773 | | CoinMpsIO m; |
774 | | m.passInMessageHandler(handler_); |
775 | | *m.messagesPointer() = coinMessages(); |
776 | | bool savePrefix = m.messageHandler()->prefix(); |
777 | | m.messageHandler()->setPrefix(handler_->prefix()); |
778 | | status = m.readBasis(fileName, "", columnActivity_, status_ + numberColumns_, |
779 | | status_, |
780 | | columnNames_, numberColumns_, |
781 | | rowNames_, numberRows_); |
782 | | m.messageHandler()->setPrefix(savePrefix); |
783 | | if (status >= 0) { |
784 | | if (!status) { |
785 | | // set values |
786 | | int iColumn, iRow; |
787 | | for (iRow = 0; iRow < numberRows_; iRow++) { |
788 | | if (getRowStatus(iRow) == atLowerBound) |
789 | | rowActivity_[iRow] = rowLower_[iRow]; |
790 | | else if (getRowStatus(iRow) == atUpperBound) |
791 | | rowActivity_[iRow] = rowUpper_[iRow]; |
792 | | } |
793 | | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
794 | | if (getColumnStatus(iColumn) == atLowerBound) |
795 | | columnActivity_[iColumn] = columnLower_[iColumn]; |
796 | | else if (getColumnStatus(iColumn) == atUpperBound) |
797 | | columnActivity_[iColumn] = columnUpper_[iColumn]; |
798 | | } |
799 | | } else { |
800 | | memset(rowActivity_, 0, numberRows_ * sizeof(double)); |
801 | | matrix_->times(-1.0, columnActivity_, rowActivity_); |
802 | | } |
803 | | } else { |
804 | | // errors |
805 | | handler_->message(CLP_IMPORT_ERRORS, messages_) |
806 | | << status << fileName << CoinMessageEol; |
807 | | } |
808 | | return status; |
| 750 | int status = 0; |
| 751 | if (strcmp(fileName, "-") != 0 && strcmp(fileName, "stdin") != 0) { |
| 752 | FILE *fp = fopen(fileName, "r"); |
| 753 | if (fp) { |
| 754 | // can open - lets go for it |
| 755 | fclose(fp); |
| 756 | } else { |
| 757 | handler_->message(CLP_UNABLE_OPEN, messages_) |
| 758 | << fileName << CoinMessageEol; |
| 759 | return -1; |
| 760 | } |
| 761 | } |
| 762 | CoinMpsIO m; |
| 763 | m.passInMessageHandler(handler_); |
| 764 | *m.messagesPointer() = coinMessages(); |
| 765 | bool savePrefix = m.messageHandler()->prefix(); |
| 766 | m.messageHandler()->setPrefix(handler_->prefix()); |
| 767 | status = m.readBasis(fileName, "", columnActivity_, status_ + numberColumns_, |
| 768 | status_, |
| 769 | columnNames_, numberColumns_, |
| 770 | rowNames_, numberRows_); |
| 771 | m.messageHandler()->setPrefix(savePrefix); |
| 772 | if (status >= 0) { |
| 773 | if (!status) { |
| 774 | // set values |
| 775 | int iColumn, iRow; |
| 776 | for (iRow = 0; iRow < numberRows_; iRow++) { |
| 777 | if (getRowStatus(iRow) == atLowerBound) |
| 778 | rowActivity_[iRow] = rowLower_[iRow]; |
| 779 | else if (getRowStatus(iRow) == atUpperBound) |
| 780 | rowActivity_[iRow] = rowUpper_[iRow]; |
| 781 | } |
| 782 | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 783 | if (getColumnStatus(iColumn) == atLowerBound) |
| 784 | columnActivity_[iColumn] = columnLower_[iColumn]; |
| 785 | else if (getColumnStatus(iColumn) == atUpperBound) |
| 786 | columnActivity_[iColumn] = columnUpper_[iColumn]; |
| 787 | } |
| 788 | } else { |
| 789 | memset(rowActivity_, 0, numberRows_ * sizeof(double)); |
| 790 | matrix_->times(-1.0, columnActivity_, rowActivity_); |
| 791 | } |
| 792 | } else { |
| 793 | // errors |
| 794 | handler_->message(CLP_IMPORT_ERRORS, messages_) |
| 795 | << status << fileName << CoinMessageEol; |
| 796 | } |
| 797 | return status; |
818 | | const ClpSimplex * model2 = static_cast<const ClpSimplex *> (this); |
819 | | bool changed = false; |
820 | | int numberChanged = 0; |
821 | | int numberFreeColumnsInPrimal=0; |
822 | | int iColumn; |
823 | | // check if we need to change bounds to rows |
824 | | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
825 | | if (columnUpper_[iColumn] < 1.0e20) { |
826 | | if (columnLower_[iColumn] > -1.0e20) { |
827 | | changed = true; |
828 | | numberChanged++; |
829 | | } |
830 | | } else if (columnLower_[iColumn] < -1.0e20) { |
831 | | numberFreeColumnsInPrimal++; |
832 | | } |
833 | | } |
834 | | int iRow; |
835 | | int numberExtraRows = 0; |
836 | | int numberFreeColumnsInDual=0; |
837 | | if (numberChanged <= fractionColumnRanges * numberColumns_) { |
838 | | for (iRow = 0; iRow < numberRows_; iRow++) { |
839 | | if (rowLower_[iRow] > -1.0e20 && |
840 | | rowUpper_[iRow] < 1.0e20) { |
841 | | if (rowUpper_[iRow] != rowLower_[iRow]) |
842 | | numberExtraRows++; |
843 | | else |
844 | | numberFreeColumnsInDual++; |
845 | | } |
846 | | } |
847 | | if (numberExtraRows > fractionRowRanges * numberRows_) |
848 | | return NULL; |
849 | | } else { |
850 | | return NULL; |
851 | | } |
852 | | printf("would have %d free columns in primal, %d in dual\n", |
853 | | numberFreeColumnsInPrimal,numberFreeColumnsInDual); |
854 | | if (4*(numberFreeColumnsInDual-numberFreeColumnsInPrimal)> |
855 | | numberColumns_&&fractionRowRanges<1.0) |
856 | | return NULL; //dangerous (well anyway in dual) |
857 | | if (changed) { |
858 | | ClpSimplex * model3 = new ClpSimplex(*model2); |
859 | | CoinBuild build; |
860 | | double one = 1.0; |
861 | | int numberColumns = model3->numberColumns(); |
862 | | const double * columnLower = model3->columnLower(); |
863 | | const double * columnUpper = model3->columnUpper(); |
864 | | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
865 | | if (columnUpper[iColumn] < 1.0e20 && |
866 | | columnLower[iColumn] > -1.0e20) { |
867 | | if (fabs(columnLower[iColumn]) < fabs(columnUpper[iColumn])) { |
868 | | double value = columnUpper[iColumn]; |
869 | | model3->setColumnUpper(iColumn, COIN_DBL_MAX); |
870 | | build.addRow(1, &iColumn, &one, -COIN_DBL_MAX, value); |
871 | | } else { |
872 | | double value = columnLower[iColumn]; |
873 | | model3->setColumnLower(iColumn, -COIN_DBL_MAX); |
874 | | build.addRow(1, &iColumn, &one, value, COIN_DBL_MAX); |
875 | | } |
876 | | } |
877 | | } |
878 | | model3->addRows(build); |
879 | | model2 = model3; |
880 | | } |
881 | | int numberColumns = model2->numberColumns(); |
882 | | const double * columnLower = model2->columnLower(); |
883 | | const double * columnUpper = model2->columnUpper(); |
884 | | int numberRows = model2->numberRows(); |
885 | | double * rowLower = CoinCopyOfArray(model2->rowLower(), numberRows); |
886 | | double * rowUpper = CoinCopyOfArray(model2->rowUpper(), numberRows); |
| 807 | const ClpSimplex *model2 = static_cast< const ClpSimplex * >(this); |
| 808 | bool changed = false; |
| 809 | int numberChanged = 0; |
| 810 | int numberFreeColumnsInPrimal = 0; |
| 811 | int iColumn; |
| 812 | // check if we need to change bounds to rows |
| 813 | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 814 | if (columnUpper_[iColumn] < 1.0e20) { |
| 815 | if (columnLower_[iColumn] > -1.0e20) { |
| 816 | changed = true; |
| 817 | numberChanged++; |
| 818 | } |
| 819 | } else if (columnLower_[iColumn] < -1.0e20) { |
| 820 | numberFreeColumnsInPrimal++; |
| 821 | } |
| 822 | } |
| 823 | int iRow; |
| 824 | int numberExtraRows = 0; |
| 825 | int numberFreeColumnsInDual = 0; |
| 826 | if (numberChanged <= fractionColumnRanges * numberColumns_) { |
| 827 | for (iRow = 0; iRow < numberRows_; iRow++) { |
| 828 | if (rowLower_[iRow] > -1.0e20 && rowUpper_[iRow] < 1.0e20) { |
| 829 | if (rowUpper_[iRow] != rowLower_[iRow]) |
| 830 | numberExtraRows++; |
| 831 | else |
| 832 | numberFreeColumnsInDual++; |
| 833 | } |
| 834 | } |
| 835 | if (numberExtraRows > fractionRowRanges * numberRows_) |
| 836 | return NULL; |
| 837 | } else { |
| 838 | return NULL; |
| 839 | } |
| 840 | printf("would have %d free columns in primal, %d in dual\n", |
| 841 | numberFreeColumnsInPrimal, numberFreeColumnsInDual); |
| 842 | if (4 * (numberFreeColumnsInDual - numberFreeColumnsInPrimal) > numberColumns_ && fractionRowRanges < 1.0) |
| 843 | return NULL; //dangerous (well anyway in dual) |
| 844 | if (changed) { |
| 845 | ClpSimplex *model3 = new ClpSimplex(*model2); |
| 846 | CoinBuild build; |
| 847 | double one = 1.0; |
| 848 | int numberColumns = model3->numberColumns(); |
| 849 | const double *columnLower = model3->columnLower(); |
| 850 | const double *columnUpper = model3->columnUpper(); |
| 851 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 852 | if (columnUpper[iColumn] < 1.0e20 && columnLower[iColumn] > -1.0e20) { |
| 853 | if (fabs(columnLower[iColumn]) < fabs(columnUpper[iColumn])) { |
| 854 | double value = columnUpper[iColumn]; |
| 855 | model3->setColumnUpper(iColumn, COIN_DBL_MAX); |
| 856 | build.addRow(1, &iColumn, &one, -COIN_DBL_MAX, value); |
| 857 | } else { |
| 858 | double value = columnLower[iColumn]; |
| 859 | model3->setColumnLower(iColumn, -COIN_DBL_MAX); |
| 860 | build.addRow(1, &iColumn, &one, value, COIN_DBL_MAX); |
| 861 | } |
| 862 | } |
| 863 | } |
| 864 | model3->addRows(build); |
| 865 | model2 = model3; |
| 866 | } |
| 867 | int numberColumns = model2->numberColumns(); |
| 868 | const double *columnLower = model2->columnLower(); |
| 869 | const double *columnUpper = model2->columnUpper(); |
| 870 | int numberRows = model2->numberRows(); |
| 871 | double *rowLower = CoinCopyOfArray(model2->rowLower(), numberRows); |
| 872 | double *rowUpper = CoinCopyOfArray(model2->rowUpper(), numberRows); |
888 | | const double * objective = model2->objective(); |
889 | | CoinPackedMatrix * matrix = model2->matrix(); |
890 | | // get transpose |
891 | | CoinPackedMatrix rowCopy = *matrix; |
892 | | const int * row = matrix->getIndices(); |
893 | | const int * columnLength = matrix->getVectorLengths(); |
894 | | const CoinBigIndex * columnStart = matrix->getVectorStarts(); |
895 | | const double * elementByColumn = matrix->getElements(); |
896 | | double objOffset = 0.0; |
897 | | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
898 | | double offset = 0.0; |
899 | | double objValue = optimizationDirection_ * objective[iColumn]; |
900 | | if (columnUpper[iColumn] > 1.0e20) { |
901 | | if (columnLower[iColumn] > -1.0e20) |
902 | | offset = columnLower[iColumn]; |
903 | | } else if (columnLower[iColumn] < -1.0e20) { |
904 | | offset = columnUpper[iColumn]; |
905 | | } else { |
906 | | // taken care of before |
907 | | abort(); |
908 | | } |
909 | | if (offset) { |
910 | | objOffset += offset * objValue; |
911 | | for (CoinBigIndex j = columnStart[iColumn]; |
912 | | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
913 | | int iRow = row[j]; |
914 | | if (rowLower[iRow] > -1.0e20) |
915 | | rowLower[iRow] -= offset * elementByColumn[j]; |
916 | | if (rowUpper[iRow] < 1.0e20) |
917 | | rowUpper[iRow] -= offset * elementByColumn[j]; |
918 | | } |
919 | | } |
920 | | } |
921 | | int * which = new int[numberRows+numberExtraRows]; |
922 | | rowCopy.reverseOrdering(); |
923 | | rowCopy.transpose(); |
924 | | double * fromRowsLower = new double[numberRows+numberExtraRows]; |
925 | | double * fromRowsUpper = new double[numberRows+numberExtraRows]; |
926 | | double * newObjective = new double[numberRows+numberExtraRows]; |
927 | | double * fromColumnsLower = new double[numberColumns]; |
928 | | double * fromColumnsUpper = new double[numberColumns]; |
929 | | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
930 | | double objValue = optimizationDirection_ * objective[iColumn]; |
931 | | // Offset is already in |
932 | | if (columnUpper[iColumn] > 1.0e20) { |
933 | | if (columnLower[iColumn] > -1.0e20) { |
934 | | fromColumnsLower[iColumn] = -COIN_DBL_MAX; |
935 | | fromColumnsUpper[iColumn] = objValue; |
936 | | } else { |
937 | | // free |
938 | | fromColumnsLower[iColumn] = objValue; |
939 | | fromColumnsUpper[iColumn] = objValue; |
940 | | } |
941 | | } else if (columnLower[iColumn] < -1.0e20) { |
942 | | fromColumnsLower[iColumn] = objValue; |
943 | | fromColumnsUpper[iColumn] = COIN_DBL_MAX; |
944 | | } else { |
945 | | abort(); |
946 | | } |
947 | | } |
948 | | int kRow = 0; |
949 | | int kExtraRow = numberRows; |
950 | | for (iRow = 0; iRow < numberRows; iRow++) { |
951 | | if (rowLower[iRow] < -1.0e20) { |
952 | | assert (rowUpper[iRow] < 1.0e20); |
953 | | newObjective[kRow] = -rowUpper[iRow]; |
954 | | fromRowsLower[kRow] = -COIN_DBL_MAX; |
955 | | fromRowsUpper[kRow] = 0.0; |
956 | | which[kRow] = iRow; |
957 | | kRow++; |
958 | | } else if (rowUpper[iRow] > 1.0e20) { |
959 | | newObjective[kRow] = -rowLower[iRow]; |
960 | | fromRowsLower[kRow] = 0.0; |
961 | | fromRowsUpper[kRow] = COIN_DBL_MAX; |
962 | | which[kRow] = iRow; |
963 | | kRow++; |
964 | | } else { |
965 | | if (rowUpper[iRow] == rowLower[iRow]) { |
966 | | newObjective[kRow] = -rowLower[iRow]; |
967 | | fromRowsLower[kRow] = -COIN_DBL_MAX;; |
968 | | fromRowsUpper[kRow] = COIN_DBL_MAX; |
969 | | which[kRow] = iRow; |
970 | | kRow++; |
971 | | } else { |
972 | | // range |
973 | | newObjective[kRow] = -rowUpper[iRow]; |
974 | | fromRowsLower[kRow] = -COIN_DBL_MAX; |
975 | | fromRowsUpper[kRow] = 0.0; |
976 | | which[kRow] = iRow; |
977 | | kRow++; |
978 | | newObjective[kExtraRow] = -rowLower[iRow]; |
979 | | fromRowsLower[kExtraRow] = 0.0; |
980 | | fromRowsUpper[kExtraRow] = COIN_DBL_MAX; |
981 | | which[kExtraRow] = iRow; |
982 | | kExtraRow++; |
983 | | } |
984 | | } |
985 | | } |
986 | | if (numberExtraRows) { |
987 | | CoinPackedMatrix newCopy; |
988 | | newCopy.setExtraGap(0.0); |
989 | | newCopy.setExtraMajor(0.0); |
990 | | newCopy.submatrixOfWithDuplicates(rowCopy, kExtraRow, which); |
991 | | rowCopy = newCopy; |
992 | | } |
993 | | ClpSimplex * modelDual = new ClpSimplex(); |
994 | | modelDual->passInEventHandler(eventHandler_); |
995 | | modelDual->loadProblem(rowCopy, fromRowsLower, fromRowsUpper, newObjective, |
996 | | fromColumnsLower, fromColumnsUpper); |
997 | | modelDual->setObjectiveOffset(objOffset); |
998 | | modelDual->setDualBound(model2->dualBound()); |
999 | | modelDual->setInfeasibilityCost(model2->infeasibilityCost()); |
1000 | | modelDual->setDualTolerance(model2->dualTolerance()); |
1001 | | modelDual->setPrimalTolerance(model2->primalTolerance()); |
1002 | | modelDual->setPerturbation(model2->perturbation()); |
1003 | | modelDual->setSpecialOptions(model2->specialOptions()); |
1004 | | modelDual->setMoreSpecialOptions(model2->moreSpecialOptions()); |
1005 | | modelDual->setMaximumIterations(model2->maximumIterations()); |
1006 | | modelDual->setFactorizationFrequency(model2->factorizationFrequency()); |
1007 | | modelDual->setLogLevel(model2->logLevel()); |
1008 | | delete [] fromRowsLower; |
1009 | | delete [] fromRowsUpper; |
1010 | | delete [] fromColumnsLower; |
1011 | | delete [] fromColumnsUpper; |
1012 | | delete [] newObjective; |
1013 | | delete [] which; |
1014 | | delete [] rowLower; |
1015 | | delete [] rowUpper; |
1016 | | if (changed) |
1017 | | delete model2; |
1018 | | modelDual->createStatus(); |
1019 | | return modelDual; |
| 874 | const double *objective = model2->objective(); |
| 875 | CoinPackedMatrix *matrix = model2->matrix(); |
| 876 | // get transpose |
| 877 | CoinPackedMatrix rowCopy = *matrix; |
| 878 | const int *row = matrix->getIndices(); |
| 879 | const int *columnLength = matrix->getVectorLengths(); |
| 880 | const CoinBigIndex *columnStart = matrix->getVectorStarts(); |
| 881 | const double *elementByColumn = matrix->getElements(); |
| 882 | double objOffset = 0.0; |
| 883 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 884 | double offset = 0.0; |
| 885 | double objValue = optimizationDirection_ * objective[iColumn]; |
| 886 | if (columnUpper[iColumn] > 1.0e20) { |
| 887 | if (columnLower[iColumn] > -1.0e20) |
| 888 | offset = columnLower[iColumn]; |
| 889 | } else if (columnLower[iColumn] < -1.0e20) { |
| 890 | offset = columnUpper[iColumn]; |
| 891 | } else { |
| 892 | // taken care of before |
| 893 | abort(); |
| 894 | } |
| 895 | if (offset) { |
| 896 | objOffset += offset * objValue; |
| 897 | for (CoinBigIndex j = columnStart[iColumn]; |
| 898 | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
| 899 | int iRow = row[j]; |
| 900 | if (rowLower[iRow] > -1.0e20) |
| 901 | rowLower[iRow] -= offset * elementByColumn[j]; |
| 902 | if (rowUpper[iRow] < 1.0e20) |
| 903 | rowUpper[iRow] -= offset * elementByColumn[j]; |
| 904 | } |
| 905 | } |
| 906 | } |
| 907 | int *which = new int[numberRows + numberExtraRows]; |
| 908 | rowCopy.reverseOrdering(); |
| 909 | rowCopy.transpose(); |
| 910 | double *fromRowsLower = new double[numberRows + numberExtraRows]; |
| 911 | double *fromRowsUpper = new double[numberRows + numberExtraRows]; |
| 912 | double *newObjective = new double[numberRows + numberExtraRows]; |
| 913 | double *fromColumnsLower = new double[numberColumns]; |
| 914 | double *fromColumnsUpper = new double[numberColumns]; |
| 915 | for (iColumn = 0; iColumn < numberColumns; iColumn++) { |
| 916 | double objValue = optimizationDirection_ * objective[iColumn]; |
| 917 | // Offset is already in |
| 918 | if (columnUpper[iColumn] > 1.0e20) { |
| 919 | if (columnLower[iColumn] > -1.0e20) { |
| 920 | fromColumnsLower[iColumn] = -COIN_DBL_MAX; |
| 921 | fromColumnsUpper[iColumn] = objValue; |
| 922 | } else { |
| 923 | // free |
| 924 | fromColumnsLower[iColumn] = objValue; |
| 925 | fromColumnsUpper[iColumn] = objValue; |
| 926 | } |
| 927 | } else if (columnLower[iColumn] < -1.0e20) { |
| 928 | fromColumnsLower[iColumn] = objValue; |
| 929 | fromColumnsUpper[iColumn] = COIN_DBL_MAX; |
| 930 | } else { |
| 931 | abort(); |
| 932 | } |
| 933 | } |
| 934 | int kRow = 0; |
| 935 | int kExtraRow = numberRows; |
| 936 | for (iRow = 0; iRow < numberRows; iRow++) { |
| 937 | if (rowLower[iRow] < -1.0e20) { |
| 938 | assert(rowUpper[iRow] < 1.0e20); |
| 939 | newObjective[kRow] = -rowUpper[iRow]; |
| 940 | fromRowsLower[kRow] = -COIN_DBL_MAX; |
| 941 | fromRowsUpper[kRow] = 0.0; |
| 942 | which[kRow] = iRow; |
| 943 | kRow++; |
| 944 | } else if (rowUpper[iRow] > 1.0e20) { |
| 945 | newObjective[kRow] = -rowLower[iRow]; |
| 946 | fromRowsLower[kRow] = 0.0; |
| 947 | fromRowsUpper[kRow] = COIN_DBL_MAX; |
| 948 | which[kRow] = iRow; |
| 949 | kRow++; |
| 950 | } else { |
| 951 | if (rowUpper[iRow] == rowLower[iRow]) { |
| 952 | newObjective[kRow] = -rowLower[iRow]; |
| 953 | fromRowsLower[kRow] = -COIN_DBL_MAX; |
| 954 | ; |
| 955 | fromRowsUpper[kRow] = COIN_DBL_MAX; |
| 956 | which[kRow] = iRow; |
| 957 | kRow++; |
| 958 | } else { |
| 959 | // range |
| 960 | newObjective[kRow] = -rowUpper[iRow]; |
| 961 | fromRowsLower[kRow] = -COIN_DBL_MAX; |
| 962 | fromRowsUpper[kRow] = 0.0; |
| 963 | which[kRow] = iRow; |
| 964 | kRow++; |
| 965 | newObjective[kExtraRow] = -rowLower[iRow]; |
| 966 | fromRowsLower[kExtraRow] = 0.0; |
| 967 | fromRowsUpper[kExtraRow] = COIN_DBL_MAX; |
| 968 | which[kExtraRow] = iRow; |
| 969 | kExtraRow++; |
| 970 | } |
| 971 | } |
| 972 | } |
| 973 | if (numberExtraRows) { |
| 974 | CoinPackedMatrix newCopy; |
| 975 | newCopy.setExtraGap(0.0); |
| 976 | newCopy.setExtraMajor(0.0); |
| 977 | newCopy.submatrixOfWithDuplicates(rowCopy, kExtraRow, which); |
| 978 | rowCopy = newCopy; |
| 979 | } |
| 980 | ClpSimplex *modelDual = new ClpSimplex(); |
| 981 | modelDual->passInEventHandler(eventHandler_); |
| 982 | modelDual->loadProblem(rowCopy, fromRowsLower, fromRowsUpper, newObjective, |
| 983 | fromColumnsLower, fromColumnsUpper); |
| 984 | modelDual->setObjectiveOffset(objOffset); |
| 985 | modelDual->setDualBound(model2->dualBound()); |
| 986 | modelDual->setInfeasibilityCost(model2->infeasibilityCost()); |
| 987 | modelDual->setDualTolerance(model2->dualTolerance()); |
| 988 | modelDual->setPrimalTolerance(model2->primalTolerance()); |
| 989 | modelDual->setPerturbation(model2->perturbation()); |
| 990 | modelDual->setSpecialOptions(model2->specialOptions()); |
| 991 | modelDual->setMoreSpecialOptions(model2->moreSpecialOptions()); |
| 992 | modelDual->setMaximumIterations(model2->maximumIterations()); |
| 993 | modelDual->setFactorizationFrequency(model2->factorizationFrequency()); |
| 994 | modelDual->setLogLevel(model2->logLevel()); |
| 995 | delete[] fromRowsLower; |
| 996 | delete[] fromRowsUpper; |
| 997 | delete[] fromColumnsLower; |
| 998 | delete[] fromColumnsUpper; |
| 999 | delete[] newObjective; |
| 1000 | delete[] which; |
| 1001 | delete[] rowLower; |
| 1002 | delete[] rowUpper; |
| 1003 | if (changed) |
| 1004 | delete model2; |
| 1005 | modelDual->createStatus(); |
| 1006 | return modelDual; |
1077 | | // position at bound information |
1078 | | int jColumn = numberRows_; |
1079 | | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
1080 | | double objValue = optimizationDirection_ * objective[iColumn]; |
1081 | | Status status = dualProblem->getRowStatus(iColumn); |
1082 | | double otherValue = COIN_DBL_MAX; |
1083 | | if (columnUpper_[iColumn] < 1.0e20 && |
1084 | | columnLower_[iColumn] > -1.0e20) { |
1085 | | if (fabs(columnLower_[iColumn]) < fabs(columnUpper_[iColumn])) { |
1086 | | otherValue = columnUpper_[iColumn] + dualDj[jColumn]; |
1087 | | } else { |
1088 | | otherValue = columnLower_[iColumn] + dualDj[jColumn]; |
1089 | | } |
1090 | | jColumn++; |
1091 | | } |
1092 | | if (status == basic) { |
1093 | | // column is at bound |
1094 | | if (otherValue == COIN_DBL_MAX) { |
1095 | | reducedCost_[iColumn] = objValue - dualActs[iColumn]; |
1096 | | if (columnUpper_[iColumn] > 1.0e20) { |
1097 | | if (columnLower_[iColumn] > -1.0e20) { |
1098 | | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
1099 | | setColumnStatus(iColumn, atLowerBound); |
1100 | | else |
1101 | | setColumnStatus(iColumn, isFixed); |
1102 | | columnActivity_[iColumn] = columnLower_[iColumn]; |
1103 | | } else { |
1104 | | // free |
1105 | | setColumnStatus(iColumn, isFree); |
1106 | | columnActivity_[iColumn] = 0.0; |
1107 | | } |
1108 | | } else { |
1109 | | setColumnStatus(iColumn, atUpperBound); |
1110 | | columnActivity_[iColumn] = columnUpper_[iColumn]; |
1111 | | } |
1112 | | } else { |
1113 | | reducedCost_[iColumn] = objValue - dualActs[iColumn]; |
1114 | | //printf("other dual sol %g\n",otherValue); |
1115 | | if (fabs(otherValue - columnLower_[iColumn]) < 1.0e-5) { |
1116 | | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
1117 | | setColumnStatus(iColumn, atLowerBound); |
1118 | | else |
1119 | | setColumnStatus(iColumn, isFixed); |
1120 | | columnActivity_[iColumn] = columnLower_[iColumn]; |
1121 | | } else if (fabs(otherValue - columnUpper_[iColumn]) < 1.0e-5) { |
1122 | | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
1123 | | setColumnStatus(iColumn, atUpperBound); |
1124 | | else |
1125 | | setColumnStatus(iColumn, isFixed); |
1126 | | columnActivity_[iColumn] = columnUpper_[iColumn]; |
1127 | | } else { |
1128 | | setColumnStatus(iColumn, superBasic); |
1129 | | columnActivity_[iColumn] = otherValue; |
1130 | | } |
1131 | | } |
| 1063 | // position at bound information |
| 1064 | int jColumn = numberRows_; |
| 1065 | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 1066 | double objValue = optimizationDirection_ * objective[iColumn]; |
| 1067 | Status status = dualProblem->getRowStatus(iColumn); |
| 1068 | double otherValue = COIN_DBL_MAX; |
| 1069 | if (columnUpper_[iColumn] < 1.0e20 && columnLower_[iColumn] > -1.0e20) { |
| 1070 | if (fabs(columnLower_[iColumn]) < fabs(columnUpper_[iColumn])) { |
| 1071 | otherValue = columnUpper_[iColumn] + dualDj[jColumn]; |
| 1072 | } else { |
| 1073 | otherValue = columnLower_[iColumn] + dualDj[jColumn]; |
| 1074 | } |
| 1075 | jColumn++; |
| 1076 | } |
| 1077 | if (status == basic) { |
| 1078 | // column is at bound |
| 1079 | if (otherValue == COIN_DBL_MAX) { |
| 1080 | reducedCost_[iColumn] = objValue - dualActs[iColumn]; |
| 1081 | if (columnUpper_[iColumn] > 1.0e20) { |
| 1082 | if (columnLower_[iColumn] > -1.0e20) { |
| 1083 | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
| 1084 | setColumnStatus(iColumn, atLowerBound); |
| 1085 | else |
| 1086 | setColumnStatus(iColumn, isFixed); |
| 1087 | columnActivity_[iColumn] = columnLower_[iColumn]; |
1133 | | if (otherValue == COIN_DBL_MAX) { |
1134 | | // column basic |
1135 | | setColumnStatus(iColumn, basic); |
1136 | | numberBasic++; |
1137 | | if (columnLower_[iColumn] > -1.0e20) { |
1138 | | columnActivity_[iColumn] = -dualDual[iColumn] + columnLower_[iColumn]; |
1139 | | } else if (columnUpper_[iColumn] < 1.0e20) { |
1140 | | columnActivity_[iColumn] = -dualDual[iColumn] + columnUpper_[iColumn]; |
1141 | | } else { |
1142 | | columnActivity_[iColumn] = -dualDual[iColumn]; |
1143 | | } |
1144 | | reducedCost_[iColumn] = 0.0; |
1145 | | } else { |
1146 | | // may be at other bound |
1147 | | //printf("xx %d %g jcol %d\n",iColumn,otherValue,jColumn-1); |
1148 | | if (dualProblem->getColumnStatus(jColumn - 1) != basic) { |
1149 | | // column basic |
1150 | | setColumnStatus(iColumn, basic); |
1151 | | numberBasic++; |
1152 | | //printf("Col %d otherV %g dualDual %g\n",iColumn, |
1153 | | // otherValue,dualDual[iColumn]); |
1154 | | columnActivity_[iColumn] = -dualDual[iColumn]; |
1155 | | columnActivity_[iColumn] = otherValue; |
1156 | | reducedCost_[iColumn] = 0.0; |
1157 | | } else { |
1158 | | reducedCost_[iColumn] = objValue - dualActs[iColumn]; |
1159 | | if (fabs(otherValue - columnLower_[iColumn]) < 1.0e-5) { |
1160 | | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
1161 | | setColumnStatus(iColumn, atLowerBound); |
1162 | | else |
1163 | | setColumnStatus(iColumn, isFixed); |
1164 | | columnActivity_[iColumn] = columnLower_[iColumn]; |
1165 | | } else if (fabs(otherValue - columnUpper_[iColumn]) < 1.0e-5) { |
1166 | | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
1167 | | setColumnStatus(iColumn, atUpperBound); |
1168 | | else |
1169 | | setColumnStatus(iColumn, isFixed); |
1170 | | columnActivity_[iColumn] = columnUpper_[iColumn]; |
1171 | | } else { |
1172 | | setColumnStatus(iColumn, superBasic); |
1173 | | columnActivity_[iColumn] = otherValue; |
1174 | | } |
1175 | | } |
1176 | | } |
1177 | | } |
1178 | | } |
1179 | | // now rows |
1180 | | int kExtraRow = jColumn; |
1181 | | int numberRanges = 0; |
1182 | | for (iRow = 0; iRow < numberRows_; iRow++) { |
1183 | | Status status = dualProblem->getColumnStatus(iRow); |
1184 | | if (status == basic) { |
1185 | | // row is at bound |
1186 | | dual_[iRow] = dualSol[iRow];; |
| 1089 | // free |
| 1090 | setColumnStatus(iColumn, isFree); |
| 1091 | columnActivity_[iColumn] = 0.0; |
| 1092 | } |
| 1093 | } else { |
| 1094 | setColumnStatus(iColumn, atUpperBound); |
| 1095 | columnActivity_[iColumn] = columnUpper_[iColumn]; |
| 1096 | } |
| 1097 | } else { |
| 1098 | reducedCost_[iColumn] = objValue - dualActs[iColumn]; |
| 1099 | //printf("other dual sol %g\n",otherValue); |
| 1100 | if (fabs(otherValue - columnLower_[iColumn]) < 1.0e-5) { |
| 1101 | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
| 1102 | setColumnStatus(iColumn, atLowerBound); |
| 1103 | else |
| 1104 | setColumnStatus(iColumn, isFixed); |
| 1105 | columnActivity_[iColumn] = columnLower_[iColumn]; |
| 1106 | } else if (fabs(otherValue - columnUpper_[iColumn]) < 1.0e-5) { |
| 1107 | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
| 1108 | setColumnStatus(iColumn, atUpperBound); |
| 1109 | else |
| 1110 | setColumnStatus(iColumn, isFixed); |
| 1111 | columnActivity_[iColumn] = columnUpper_[iColumn]; |
| 1112 | } else { |
| 1113 | setColumnStatus(iColumn, superBasic); |
| 1114 | columnActivity_[iColumn] = otherValue; |
| 1115 | } |
| 1116 | } |
| 1117 | } else { |
| 1118 | if (otherValue == COIN_DBL_MAX) { |
| 1119 | // column basic |
| 1120 | setColumnStatus(iColumn, basic); |
| 1121 | numberBasic++; |
| 1122 | if (columnLower_[iColumn] > -1.0e20) { |
| 1123 | columnActivity_[iColumn] = -dualDual[iColumn] + columnLower_[iColumn]; |
| 1124 | } else if (columnUpper_[iColumn] < 1.0e20) { |
| 1125 | columnActivity_[iColumn] = -dualDual[iColumn] + columnUpper_[iColumn]; |
| 1126 | } else { |
| 1127 | columnActivity_[iColumn] = -dualDual[iColumn]; |
| 1128 | } |
| 1129 | reducedCost_[iColumn] = 0.0; |
| 1130 | } else { |
| 1131 | // may be at other bound |
| 1132 | //printf("xx %d %g jcol %d\n",iColumn,otherValue,jColumn-1); |
| 1133 | if (dualProblem->getColumnStatus(jColumn - 1) != basic) { |
| 1134 | // column basic |
| 1135 | setColumnStatus(iColumn, basic); |
| 1136 | numberBasic++; |
| 1137 | //printf("Col %d otherV %g dualDual %g\n",iColumn, |
| 1138 | // otherValue,dualDual[iColumn]); |
| 1139 | columnActivity_[iColumn] = -dualDual[iColumn]; |
| 1140 | columnActivity_[iColumn] = otherValue; |
| 1141 | reducedCost_[iColumn] = 0.0; |
| 1142 | } else { |
| 1143 | reducedCost_[iColumn] = objValue - dualActs[iColumn]; |
| 1144 | if (fabs(otherValue - columnLower_[iColumn]) < 1.0e-5) { |
| 1145 | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
| 1146 | setColumnStatus(iColumn, atLowerBound); |
| 1147 | else |
| 1148 | setColumnStatus(iColumn, isFixed); |
| 1149 | columnActivity_[iColumn] = columnLower_[iColumn]; |
| 1150 | } else if (fabs(otherValue - columnUpper_[iColumn]) < 1.0e-5) { |
| 1151 | if (columnUpper_[iColumn] > columnLower_[iColumn]) |
| 1152 | setColumnStatus(iColumn, atUpperBound); |
| 1153 | else |
| 1154 | setColumnStatus(iColumn, isFixed); |
| 1155 | columnActivity_[iColumn] = columnUpper_[iColumn]; |
1188 | | // row basic |
1189 | | setRowStatus(iRow, basic); |
1190 | | numberBasic++; |
1191 | | dual_[iRow] = 0.0; |
1192 | | } |
1193 | | if (rowLower_[iRow] < -1.0e20) { |
1194 | | if (status == basic) { |
1195 | | rowActivity_[iRow] = rowUpper_[iRow]; |
1196 | | setRowStatus(iRow, atUpperBound); |
1197 | | } else { |
1198 | | // might be stopped assert (dualDj[iRow] < 1.0e-5); |
1199 | | rowActivity_[iRow] = rowUpper_[iRow] + dualDj[iRow]; |
1200 | | } |
1201 | | } else if (rowUpper_[iRow] > 1.0e20) { |
1202 | | if (status == basic) { |
1203 | | rowActivity_[iRow] = rowLower_[iRow]; |
1204 | | setRowStatus(iRow, atLowerBound); |
1205 | | } else { |
1206 | | rowActivity_[iRow] = rowLower_[iRow] + dualDj[iRow]; |
1207 | | // might be stopped assert (dualDj[iRow] > -1.0e-5); |
1208 | | } |
1209 | | } else { |
1210 | | if (rowUpper_[iRow] == rowLower_[iRow]) { |
1211 | | rowActivity_[iRow] = rowLower_[iRow]; |
1212 | | if (status == basic) { |
1213 | | setRowStatus(iRow, isFixed); |
1214 | | } |
1215 | | } else { |
1216 | | // range |
1217 | | numberRanges++; |
1218 | | Status statusL = dualProblem->getColumnStatus(kExtraRow); |
1219 | | //printf("range row %d (%d), extra %d (%d) - dualSol %g,%g dualDj %g,%g\n", |
1220 | | // iRow,status,kExtraRow,statusL, dualSol[iRow], |
1221 | | // dualSol[kExtraRow],dualDj[iRow],dualDj[kExtraRow]); |
1222 | | if (status == basic) { |
1223 | | // might be stopped assert (statusL != basic); |
1224 | | rowActivity_[iRow] = rowUpper_[iRow]; |
1225 | | setRowStatus(iRow, atUpperBound); |
1226 | | } else if (statusL == basic) { |
1227 | | numberBasic--; // already counted |
1228 | | rowActivity_[iRow] = rowLower_[iRow]; |
1229 | | setRowStatus(iRow, atLowerBound); |
1230 | | dual_[iRow] = dualSol[kExtraRow];; |
1231 | | } else { |
1232 | | rowActivity_[iRow] = rowLower_[iRow] - dualDj[iRow]; |
1233 | | // might be stopped assert (dualDj[iRow] < 1.0e-5); |
1234 | | // row basic |
1235 | | //setRowStatus(iRow,basic); |
1236 | | //numberBasic++; |
1237 | | dual_[iRow] = 0.0; |
1238 | | } |
1239 | | kExtraRow++; |
1240 | | } |
1241 | | } |
1242 | | } |
1243 | | if (numberBasic != numberRows_ && 0) { |
1244 | | printf("Bad basis - ranges - coding needed\n"); |
1245 | | assert (numberRanges); |
1246 | | abort(); |
1247 | | } |
1248 | | if (optimizationDirection_ < 0.0) { |
1249 | | for (iRow = 0; iRow < numberRows_; iRow++) { |
1250 | | dual_[iRow] = -dual_[iRow]; |
1251 | | } |
1252 | | } |
1253 | | // redo row activities |
1254 | | memset(rowActivity_, 0, numberRows_ * sizeof(double)); |
1255 | | matrix_->times(1.0, columnActivity_, rowActivity_); |
1256 | | // redo reduced costs |
1257 | | memcpy(reducedCost_, this->objective(), numberColumns_ * sizeof(double)); |
1258 | | matrix_->transposeTimes(-1.0, dual_, reducedCost_); |
1259 | | checkSolutionInternal(); |
1260 | | if (sumDualInfeasibilities_ > 1.0e-5 || sumPrimalInfeasibilities_ > 1.0e-5) { |
1261 | | returnCode = 1; |
| 1157 | setColumnStatus(iColumn, superBasic); |
| 1158 | columnActivity_[iColumn] = otherValue; |
| 1159 | } |
| 1160 | } |
| 1161 | } |
| 1162 | } |
| 1163 | } |
| 1164 | // now rows |
| 1165 | int kExtraRow = jColumn; |
| 1166 | int numberRanges = 0; |
| 1167 | for (iRow = 0; iRow < numberRows_; iRow++) { |
| 1168 | Status status = dualProblem->getColumnStatus(iRow); |
| 1169 | if (status == basic) { |
| 1170 | // row is at bound |
| 1171 | dual_[iRow] = dualSol[iRow]; |
| 1172 | ; |
| 1173 | } else { |
| 1174 | // row basic |
| 1175 | setRowStatus(iRow, basic); |
| 1176 | numberBasic++; |
| 1177 | dual_[iRow] = 0.0; |
| 1178 | } |
| 1179 | if (rowLower_[iRow] < -1.0e20) { |
| 1180 | if (status == basic) { |
| 1181 | rowActivity_[iRow] = rowUpper_[iRow]; |
| 1182 | setRowStatus(iRow, atUpperBound); |
| 1183 | } else { |
| 1184 | // might be stopped assert (dualDj[iRow] < 1.0e-5); |
| 1185 | rowActivity_[iRow] = rowUpper_[iRow] + dualDj[iRow]; |
| 1186 | } |
| 1187 | } else if (rowUpper_[iRow] > 1.0e20) { |
| 1188 | if (status == basic) { |
| 1189 | rowActivity_[iRow] = rowLower_[iRow]; |
| 1190 | setRowStatus(iRow, atLowerBound); |
| 1191 | } else { |
| 1192 | rowActivity_[iRow] = rowLower_[iRow] + dualDj[iRow]; |
| 1193 | // might be stopped assert (dualDj[iRow] > -1.0e-5); |
| 1194 | } |
| 1195 | } else { |
| 1196 | if (rowUpper_[iRow] == rowLower_[iRow]) { |
| 1197 | rowActivity_[iRow] = rowLower_[iRow]; |
| 1198 | if (status == basic) { |
| 1199 | setRowStatus(iRow, isFixed); |
| 1200 | } |
| 1201 | } else { |
| 1202 | // range |
| 1203 | numberRanges++; |
| 1204 | Status statusL = dualProblem->getColumnStatus(kExtraRow); |
| 1205 | //printf("range row %d (%d), extra %d (%d) - dualSol %g,%g dualDj %g,%g\n", |
| 1206 | // iRow,status,kExtraRow,statusL, dualSol[iRow], |
| 1207 | // dualSol[kExtraRow],dualDj[iRow],dualDj[kExtraRow]); |
| 1208 | if (status == basic) { |
| 1209 | // might be stopped assert (statusL != basic); |
| 1210 | rowActivity_[iRow] = rowUpper_[iRow]; |
| 1211 | setRowStatus(iRow, atUpperBound); |
| 1212 | } else if (statusL == basic) { |
| 1213 | numberBasic--; // already counted |
| 1214 | rowActivity_[iRow] = rowLower_[iRow]; |
| 1215 | setRowStatus(iRow, atLowerBound); |
| 1216 | dual_[iRow] = dualSol[kExtraRow]; |
| 1217 | ; |
| 1218 | } else { |
| 1219 | rowActivity_[iRow] = rowLower_[iRow] - dualDj[iRow]; |
| 1220 | // might be stopped assert (dualDj[iRow] < 1.0e-5); |
| 1221 | // row basic |
| 1222 | //setRowStatus(iRow,basic); |
| 1223 | //numberBasic++; |
| 1224 | dual_[iRow] = 0.0; |
| 1225 | } |
| 1226 | kExtraRow++; |
| 1227 | } |
| 1228 | } |
| 1229 | } |
| 1230 | if (numberBasic != numberRows_ && 0) { |
| 1231 | printf("Bad basis - ranges - coding needed\n"); |
| 1232 | assert(numberRanges); |
| 1233 | abort(); |
| 1234 | } |
| 1235 | if (optimizationDirection_ < 0.0) { |
| 1236 | for (iRow = 0; iRow < numberRows_; iRow++) { |
| 1237 | dual_[iRow] = -dual_[iRow]; |
| 1238 | } |
| 1239 | } |
| 1240 | // redo row activities |
| 1241 | memset(rowActivity_, 0, numberRows_ * sizeof(double)); |
| 1242 | matrix_->times(1.0, columnActivity_, rowActivity_); |
| 1243 | // redo reduced costs |
| 1244 | memcpy(reducedCost_, this->objective(), numberColumns_ * sizeof(double)); |
| 1245 | matrix_->transposeTimes(-1.0, dual_, reducedCost_); |
| 1246 | checkSolutionInternal(); |
| 1247 | if (sumDualInfeasibilities_ > 1.0e-5 || sumPrimalInfeasibilities_ > 1.0e-5) { |
| 1248 | returnCode = 1; |
1424 | | CoinZeroN(rhs, numberRows_); |
1425 | | int iColumn; |
1426 | | int iRow; |
1427 | | CoinZeroN(whichRow, numberRows_); |
1428 | | int * backColumn = whichColumn + numberColumns_; |
1429 | | int numberRows2 = 0; |
1430 | | int numberColumns2 = 0; |
1431 | | double offset = 0.0; |
1432 | | const double * objective = this->objective(); |
1433 | | double * solution = columnActivity_; |
1434 | | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
1435 | | double lower = columnLower_[iColumn]; |
1436 | | double upper = columnUpper_[iColumn]; |
1437 | | if (upper > lower || getColumnStatus(iColumn) == ClpSimplex::basic) { |
1438 | | backColumn[iColumn] = numberColumns2; |
1439 | | whichColumn[numberColumns2++] = iColumn; |
1440 | | for (CoinBigIndex j = columnStart[iColumn]; |
1441 | | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
1442 | | int iRow = row[j]; |
1443 | | int n = whichRow[iRow]; |
1444 | | if (n == 0 && element[j]) |
1445 | | whichRow[iRow] = -iColumn - 1; |
1446 | | else if (n < 0) |
1447 | | whichRow[iRow] = 2; |
1448 | | } |
1449 | | } else { |
1450 | | // fixed |
1451 | | backColumn[iColumn] = -1; |
1452 | | solution[iColumn] = upper; |
1453 | | if (upper) { |
1454 | | offset += objective[iColumn] * upper; |
1455 | | for (CoinBigIndex j = columnStart[iColumn]; |
1456 | | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
1457 | | int iRow = row[j]; |
1458 | | double value = element[j]; |
1459 | | rhs[iRow] += upper * value; |
1460 | | } |
1461 | | } |
1462 | | } |
1463 | | } |
1464 | | int returnCode = 0; |
1465 | | double tolerance = primalTolerance(); |
1466 | | nBound = 2 * numberRows_; |
1467 | | for (iRow = 0; iRow < numberRows_; iRow++) { |
1468 | | int n = whichRow[iRow]; |
1469 | | if (n > 0) { |
1470 | | whichRow[numberRows2++] = iRow; |
1471 | | } else if (n < 0) { |
1472 | | //whichRow[numberRows2++]=iRow; |
1473 | | //continue; |
1474 | | // Can only do in certain circumstances as we don't know current value |
1475 | | if (rowLower_[iRow] == rowUpper_[iRow] || getRowStatus(iRow) == ClpSimplex::basic) { |
1476 | | // save row and column for bound |
1477 | | whichRow[--nBound] = iRow; |
1478 | | whichRow[nBound+numberRows_] = -n - 1; |
1479 | | } else if (moreBounds) { |
1480 | | // save row and column for bound |
1481 | | whichRow[--nBound] = iRow; |
1482 | | whichRow[nBound+numberRows_] = -n - 1; |
1483 | | } else { |
1484 | | whichRow[numberRows2++] = iRow; |
1485 | | } |
1486 | | } else { |
1487 | | // empty |
1488 | | double rhsValue = rhs[iRow]; |
1489 | | if (rhsValue < rowLower_[iRow] - tolerance || rhsValue > rowUpper_[iRow] + tolerance) { |
1490 | | returnCode = 1; // infeasible |
1491 | | } |
1492 | | } |
1493 | | } |
1494 | | ClpSimplex * small = NULL; |
1495 | | if (!returnCode) { |
1496 | | //printf("CRUNCH from (%d,%d) to (%d,%d)\n", |
1497 | | // numberRows_,numberColumns_,numberRows2,numberColumns2); |
1498 | | small = new ClpSimplex(this, numberRows2, whichRow, |
1499 | | numberColumns2, whichColumn, true, false); |
| 1406 | CoinZeroN(rhs, numberRows_); |
| 1407 | int iColumn; |
| 1408 | int iRow; |
| 1409 | CoinZeroN(whichRow, numberRows_); |
| 1410 | int *backColumn = whichColumn + numberColumns_; |
| 1411 | int numberRows2 = 0; |
| 1412 | int numberColumns2 = 0; |
| 1413 | double offset = 0.0; |
| 1414 | const double *objective = this->objective(); |
| 1415 | double *solution = columnActivity_; |
| 1416 | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 1417 | double lower = columnLower_[iColumn]; |
| 1418 | double upper = columnUpper_[iColumn]; |
| 1419 | if (upper > lower || getColumnStatus(iColumn) == ClpSimplex::basic) { |
| 1420 | backColumn[iColumn] = numberColumns2; |
| 1421 | whichColumn[numberColumns2++] = iColumn; |
| 1422 | for (CoinBigIndex j = columnStart[iColumn]; |
| 1423 | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
| 1424 | int iRow = row[j]; |
| 1425 | int n = whichRow[iRow]; |
| 1426 | if (n == 0 && element[j]) |
| 1427 | whichRow[iRow] = -iColumn - 1; |
| 1428 | else if (n < 0) |
| 1429 | whichRow[iRow] = 2; |
| 1430 | } |
| 1431 | } else { |
| 1432 | // fixed |
| 1433 | backColumn[iColumn] = -1; |
| 1434 | solution[iColumn] = upper; |
| 1435 | if (upper) { |
| 1436 | offset += objective[iColumn] * upper; |
| 1437 | for (CoinBigIndex j = columnStart[iColumn]; |
| 1438 | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
| 1439 | int iRow = row[j]; |
| 1440 | double value = element[j]; |
| 1441 | rhs[iRow] += upper * value; |
| 1442 | } |
| 1443 | } |
| 1444 | } |
| 1445 | } |
| 1446 | int returnCode = 0; |
| 1447 | double tolerance = primalTolerance(); |
| 1448 | nBound = 2 * numberRows_; |
| 1449 | for (iRow = 0; iRow < numberRows_; iRow++) { |
| 1450 | int n = whichRow[iRow]; |
| 1451 | if (n > 0) { |
| 1452 | whichRow[numberRows2++] = iRow; |
| 1453 | } else if (n < 0) { |
| 1454 | //whichRow[numberRows2++]=iRow; |
| 1455 | //continue; |
| 1456 | // Can only do in certain circumstances as we don't know current value |
| 1457 | if (rowLower_[iRow] == rowUpper_[iRow] || getRowStatus(iRow) == ClpSimplex::basic) { |
| 1458 | // save row and column for bound |
| 1459 | whichRow[--nBound] = iRow; |
| 1460 | whichRow[nBound + numberRows_] = -n - 1; |
| 1461 | } else if (moreBounds) { |
| 1462 | // save row and column for bound |
| 1463 | whichRow[--nBound] = iRow; |
| 1464 | whichRow[nBound + numberRows_] = -n - 1; |
| 1465 | } else { |
| 1466 | whichRow[numberRows2++] = iRow; |
| 1467 | } |
| 1468 | } else { |
| 1469 | // empty |
| 1470 | double rhsValue = rhs[iRow]; |
| 1471 | if (rhsValue < rowLower_[iRow] - tolerance || rhsValue > rowUpper_[iRow] + tolerance) { |
| 1472 | returnCode = 1; // infeasible |
| 1473 | } |
| 1474 | } |
| 1475 | } |
| 1476 | ClpSimplex *small = NULL; |
| 1477 | if (!returnCode) { |
| 1478 | //printf("CRUNCH from (%d,%d) to (%d,%d)\n", |
| 1479 | // numberRows_,numberColumns_,numberRows2,numberColumns2); |
| 1480 | small = new ClpSimplex(this, numberRows2, whichRow, |
| 1481 | numberColumns2, whichColumn, true, false); |
1519 | | CoinBigIndex numberElements = getNumElements(); |
1520 | | CoinBigIndex numberElements2 = small->getNumElements(); |
1521 | | small->setObjectiveOffset(objectiveOffset() - offset); |
1522 | | handler_->message(CLP_CRUNCH_STATS, messages_) |
1523 | | << numberRows2 << -(numberRows_ - numberRows2) |
1524 | | << numberColumns2 << -(numberColumns_ - numberColumns2) |
1525 | | << numberElements2 << -(numberElements - numberElements2) |
1526 | | << CoinMessageEol; |
1527 | | // And set objective value to match |
1528 | | small->setObjectiveValue(this->objectiveValue()); |
1529 | | double * rowLower2 = small->rowLower(); |
1530 | | double * rowUpper2 = small->rowUpper(); |
1531 | | int jRow; |
1532 | | for (jRow = 0; jRow < numberRows2; jRow++) { |
1533 | | iRow = whichRow[jRow]; |
1534 | | if (rowLower2[jRow] > -1.0e20) |
1535 | | rowLower2[jRow] -= rhs[iRow]; |
1536 | | if (rowUpper2[jRow] < 1.0e20) |
1537 | | rowUpper2[jRow] -= rhs[iRow]; |
1538 | | } |
1539 | | // and bounds |
1540 | | double * columnLower2 = small->columnLower(); |
1541 | | double * columnUpper2 = small->columnUpper(); |
1542 | | const char * integerInformation = integerType_; |
1543 | | for (jRow = nBound; jRow < 2 * numberRows_; jRow++) { |
1544 | | iRow = whichRow[jRow]; |
1545 | | iColumn = whichRow[jRow+numberRows_]; |
1546 | | double lowerRow = rowLower_[iRow]; |
1547 | | if (lowerRow > -1.0e20) |
1548 | | lowerRow -= rhs[iRow]; |
1549 | | double upperRow = rowUpper_[iRow]; |
1550 | | if (upperRow < 1.0e20) |
1551 | | upperRow -= rhs[iRow]; |
1552 | | int jColumn = backColumn[iColumn]; |
1553 | | double lower = columnLower2[jColumn]; |
1554 | | double upper = columnUpper2[jColumn]; |
1555 | | double value = 0.0; |
1556 | | for (CoinBigIndex j = columnStart[iColumn]; |
1557 | | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
1558 | | if (iRow == row[j]) { |
1559 | | value = element[j]; |
1560 | | break; |
| 1501 | CoinBigIndex numberElements = getNumElements(); |
| 1502 | CoinBigIndex numberElements2 = small->getNumElements(); |
| 1503 | small->setObjectiveOffset(objectiveOffset() - offset); |
| 1504 | handler_->message(CLP_CRUNCH_STATS, messages_) |
| 1505 | << numberRows2 << -(numberRows_ - numberRows2) |
| 1506 | << numberColumns2 << -(numberColumns_ - numberColumns2) |
| 1507 | << numberElements2 << -(numberElements - numberElements2) |
| 1508 | << CoinMessageEol; |
| 1509 | // And set objective value to match |
| 1510 | small->setObjectiveValue(this->objectiveValue()); |
| 1511 | double *rowLower2 = small->rowLower(); |
| 1512 | double *rowUpper2 = small->rowUpper(); |
| 1513 | int jRow; |
| 1514 | for (jRow = 0; jRow < numberRows2; jRow++) { |
| 1515 | iRow = whichRow[jRow]; |
| 1516 | if (rowLower2[jRow] > -1.0e20) |
| 1517 | rowLower2[jRow] -= rhs[iRow]; |
| 1518 | if (rowUpper2[jRow] < 1.0e20) |
| 1519 | rowUpper2[jRow] -= rhs[iRow]; |
| 1520 | } |
| 1521 | // and bounds |
| 1522 | double *columnLower2 = small->columnLower(); |
| 1523 | double *columnUpper2 = small->columnUpper(); |
| 1524 | const char *integerInformation = integerType_; |
| 1525 | for (jRow = nBound; jRow < 2 * numberRows_; jRow++) { |
| 1526 | iRow = whichRow[jRow]; |
| 1527 | iColumn = whichRow[jRow + numberRows_]; |
| 1528 | double lowerRow = rowLower_[iRow]; |
| 1529 | if (lowerRow > -1.0e20) |
| 1530 | lowerRow -= rhs[iRow]; |
| 1531 | double upperRow = rowUpper_[iRow]; |
| 1532 | if (upperRow < 1.0e20) |
| 1533 | upperRow -= rhs[iRow]; |
| 1534 | int jColumn = backColumn[iColumn]; |
| 1535 | double lower = columnLower2[jColumn]; |
| 1536 | double upper = columnUpper2[jColumn]; |
| 1537 | double value = 0.0; |
| 1538 | for (CoinBigIndex j = columnStart[iColumn]; |
| 1539 | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
| 1540 | if (iRow == row[j]) { |
| 1541 | value = element[j]; |
| 1542 | break; |
| 1543 | } |
| 1544 | } |
| 1545 | assert(value); |
| 1546 | // convert rowLower and Upper to implied bounds on column |
| 1547 | double newLower = -COIN_DBL_MAX; |
| 1548 | double newUpper = COIN_DBL_MAX; |
| 1549 | if (value > 0.0) { |
| 1550 | if (lowerRow > -1.0e20) |
| 1551 | newLower = lowerRow / value; |
| 1552 | if (upperRow < 1.0e20) |
| 1553 | newUpper = upperRow / value; |
| 1554 | } else { |
| 1555 | if (upperRow < 1.0e20) |
| 1556 | newLower = upperRow / value; |
| 1557 | if (lowerRow > -1.0e20) |
| 1558 | newUpper = lowerRow / value; |
| 1559 | } |
| 1560 | if (integerInformation && integerInformation[iColumn]) { |
| 1561 | if (newLower - floor(newLower) < 10.0 * tolerance) |
| 1562 | newLower = floor(newLower); |
| 1563 | else |
| 1564 | newLower = ceil(newLower); |
| 1565 | if (ceil(newUpper) - newUpper < 10.0 * tolerance) |
| 1566 | newUpper = ceil(newUpper); |
| 1567 | else |
| 1568 | newUpper = floor(newUpper); |
| 1569 | } |
| 1570 | newLower = CoinMax(lower, newLower); |
| 1571 | newUpper = CoinMin(upper, newUpper); |
| 1572 | if (newLower > newUpper + tolerance) { |
| 1573 | //printf("XXYY inf on bound\n"); |
| 1574 | returnCode = 1; |
| 1575 | } |
| 1576 | columnLower2[jColumn] = newLower; |
| 1577 | columnUpper2[jColumn] = CoinMax(newLower, newUpper); |
| 1578 | if (getRowStatus(iRow) != ClpSimplex::basic) { |
| 1579 | if (getColumnStatus(iColumn) == ClpSimplex::basic) { |
| 1580 | if (columnLower2[jColumn] == columnUpper2[jColumn]) { |
| 1581 | // can only get here if will be fixed |
| 1582 | small->setColumnStatus(jColumn, ClpSimplex::isFixed); |
| 1583 | } else { |
| 1584 | // solution is valid |
| 1585 | if (fabs(columnActivity_[iColumn] - columnLower2[jColumn]) < fabs(columnActivity_[iColumn] - columnUpper2[jColumn])) |
| 1586 | small->setColumnStatus(jColumn, ClpSimplex::atLowerBound); |
| 1587 | else |
| 1588 | small->setColumnStatus(jColumn, ClpSimplex::atUpperBound); |
| 1589 | } |
| 1590 | } else { |
| 1591 | //printf("what now neither basic\n"); |
| 1592 | } |
| 1593 | } |
| 1594 | } |
| 1595 | if (returnCode) { |
| 1596 | delete small; |
| 1597 | small = NULL; |
| 1598 | } else if (tightenBounds && integerInformation) { |
| 1599 | // See if we can tighten any bounds |
| 1600 | // use rhs for upper and small duals for lower |
| 1601 | double *up = rhs; |
| 1602 | double *lo = small->dualRowSolution(); |
| 1603 | const double *element = small->clpMatrix()->getElements(); |
| 1604 | const int *row = small->clpMatrix()->getIndices(); |
| 1605 | const CoinBigIndex *columnStart = small->clpMatrix()->getVectorStarts(); |
| 1606 | //const int * columnLength = small->clpMatrix()->getVectorLengths(); |
| 1607 | CoinZeroN(lo, numberRows2); |
| 1608 | CoinZeroN(up, numberRows2); |
| 1609 | for (int iColumn = 0; iColumn < numberColumns2; iColumn++) { |
| 1610 | double upper = columnUpper2[iColumn]; |
| 1611 | double lower = columnLower2[iColumn]; |
| 1612 | //assert (columnLength[iColumn]==columnStart[iColumn+1]-columnStart[iColumn]); |
| 1613 | for (CoinBigIndex j = columnStart[iColumn]; j < columnStart[iColumn + 1]; j++) { |
| 1614 | int iRow = row[j]; |
| 1615 | double value = element[j]; |
| 1616 | if (value > 0.0) { |
| 1617 | if (upper < 1.0e20) |
| 1618 | up[iRow] += upper * value; |
| 1619 | else |
| 1620 | up[iRow] = COIN_DBL_MAX; |
| 1621 | if (lower > -1.0e20) |
| 1622 | lo[iRow] += lower * value; |
| 1623 | else |
| 1624 | lo[iRow] = -COIN_DBL_MAX; |
| 1625 | } else { |
| 1626 | if (upper < 1.0e20) |
| 1627 | lo[iRow] += upper * value; |
| 1628 | else |
| 1629 | lo[iRow] = -COIN_DBL_MAX; |
| 1630 | if (lower > -1.0e20) |
| 1631 | up[iRow] += lower * value; |
| 1632 | else |
| 1633 | up[iRow] = COIN_DBL_MAX; |
| 1634 | } |
| 1635 | } |
| 1636 | } |
| 1637 | double *rowLower2 = small->rowLower(); |
| 1638 | double *rowUpper2 = small->rowUpper(); |
| 1639 | bool feasible = true; |
| 1640 | // make safer |
| 1641 | for (int iRow = 0; iRow < numberRows2; iRow++) { |
| 1642 | double lower = lo[iRow]; |
| 1643 | if (lower > rowUpper2[iRow] + tolerance) { |
| 1644 | feasible = false; |
| 1645 | break; |
| 1646 | } else { |
| 1647 | lo[iRow] = CoinMin(lower - rowUpper2[iRow], 0.0) - tolerance; |
| 1648 | } |
| 1649 | double upper = up[iRow]; |
| 1650 | if (upper < rowLower2[iRow] - tolerance) { |
| 1651 | feasible = false; |
| 1652 | break; |
| 1653 | } else { |
| 1654 | up[iRow] = CoinMax(upper - rowLower2[iRow], 0.0) + tolerance; |
| 1655 | } |
| 1656 | // tighten row bounds |
| 1657 | if (lower > -1.0e10) |
| 1658 | rowLower2[iRow] = CoinMax(rowLower2[iRow], |
| 1659 | lower - 1.0e-6 * (1.0 + fabs(lower))); |
| 1660 | if (upper < 1.0e10) |
| 1661 | rowUpper2[iRow] = CoinMin(rowUpper2[iRow], |
| 1662 | upper + 1.0e-6 * (1.0 + fabs(upper))); |
| 1663 | } |
| 1664 | if (!feasible) { |
| 1665 | delete small; |
| 1666 | small = NULL; |
| 1667 | } else { |
| 1668 | // and tighten |
| 1669 | for (int iColumn = 0; iColumn < numberColumns2; iColumn++) { |
| 1670 | if (integerInformation[whichColumn[iColumn]]) { |
| 1671 | double upper = columnUpper2[iColumn]; |
| 1672 | double lower = columnLower2[iColumn]; |
| 1673 | double newUpper = upper; |
| 1674 | double newLower = lower; |
| 1675 | double difference = upper - lower; |
| 1676 | if (lower > -1000.0 && upper < 1000.0) { |
| 1677 | for (CoinBigIndex j = columnStart[iColumn]; j < columnStart[iColumn + 1]; j++) { |
| 1678 | int iRow = row[j]; |
| 1679 | double value = element[j]; |
| 1680 | if (value > 0.0) { |
| 1681 | double upWithOut = up[iRow] - value * difference; |
| 1682 | if (upWithOut < 0.0) { |
| 1683 | newLower = CoinMax(newLower, lower - (upWithOut + tolerance) / value); |
| 1684 | } |
| 1685 | double lowWithOut = lo[iRow] + value * difference; |
| 1686 | if (lowWithOut > 0.0) { |
| 1687 | newUpper = CoinMin(newUpper, upper - (lowWithOut - tolerance) / value); |
| 1688 | } |
| 1689 | } else { |
| 1690 | double upWithOut = up[iRow] + value * difference; |
| 1691 | if (upWithOut < 0.0) { |
| 1692 | newUpper = CoinMin(newUpper, upper - (upWithOut + tolerance) / value); |
| 1693 | } |
| 1694 | double lowWithOut = lo[iRow] - value * difference; |
| 1695 | if (lowWithOut > 0.0) { |
| 1696 | newLower = CoinMax(newLower, lower - (lowWithOut - tolerance) / value); |
| 1697 | } |
| 1698 | } |
| 1699 | } |
| 1700 | if (newLower > lower || newUpper < upper) { |
| 1701 | if (fabs(newUpper - floor(newUpper + 0.5)) > 1.0e-6) |
| 1702 | newUpper = floor(newUpper); |
| 1703 | else |
| 1704 | newUpper = floor(newUpper + 0.5); |
| 1705 | if (fabs(newLower - ceil(newLower - 0.5)) > 1.0e-6) |
| 1706 | newLower = ceil(newLower); |
| 1707 | else |
| 1708 | newLower = ceil(newLower - 0.5); |
| 1709 | // change may be too small - check |
| 1710 | if (newLower > lower || newUpper < upper) { |
| 1711 | if (newUpper >= newLower) { |
| 1712 | // Could also tighten in this |
| 1713 | //printf("%d bounds %g %g tightened to %g %g\n", |
| 1714 | // iColumn,columnLower2[iColumn],columnUpper2[iColumn], |
| 1715 | // newLower,newUpper); |
| 1716 | #if 1 |
| 1717 | columnUpper2[iColumn] = newUpper; |
| 1718 | columnLower2[iColumn] = newLower; |
| 1719 | columnUpper_[whichColumn[iColumn]] = newUpper; |
| 1720 | columnLower_[whichColumn[iColumn]] = newLower; |
| 1721 | #endif |
| 1722 | // and adjust bounds on rows |
| 1723 | newUpper -= upper; |
| 1724 | newLower -= lower; |
| 1725 | for (CoinBigIndex j = columnStart[iColumn]; j < columnStart[iColumn + 1]; j++) { |
| 1726 | int iRow = row[j]; |
| 1727 | double value = element[j]; |
| 1728 | if (value > 0.0) { |
| 1729 | up[iRow] += newUpper * value; |
| 1730 | lo[iRow] += newLower * value; |
| 1731 | } else { |
| 1732 | lo[iRow] += newUpper * value; |
| 1733 | up[iRow] += newLower * value; |
| 1734 | } |
1562 | | } |
1563 | | assert (value); |
1564 | | // convert rowLower and Upper to implied bounds on column |
1565 | | double newLower = -COIN_DBL_MAX; |
1566 | | double newUpper = COIN_DBL_MAX; |
1567 | | if (value > 0.0) { |
1568 | | if (lowerRow > -1.0e20) |
1569 | | newLower = lowerRow / value; |
1570 | | if (upperRow < 1.0e20) |
1571 | | newUpper = upperRow / value; |
1572 | | } else { |
1573 | | if (upperRow < 1.0e20) |
1574 | | newLower = upperRow / value; |
1575 | | if (lowerRow > -1.0e20) |
1576 | | newUpper = lowerRow / value; |
1577 | | } |
1578 | | if (integerInformation && integerInformation[iColumn]) { |
1579 | | if (newLower - floor(newLower) < 10.0 * tolerance) |
1580 | | newLower = floor(newLower); |
1581 | | else |
1582 | | newLower = ceil(newLower); |
1583 | | if (ceil(newUpper) - newUpper < 10.0 * tolerance) |
1584 | | newUpper = ceil(newUpper); |
1585 | | else |
1586 | | newUpper = floor(newUpper); |
1587 | | } |
1588 | | newLower = CoinMax(lower, newLower); |
1589 | | newUpper = CoinMin(upper, newUpper); |
1590 | | if (newLower > newUpper + tolerance) { |
1591 | | //printf("XXYY inf on bound\n"); |
1592 | | returnCode = 1; |
1593 | | } |
1594 | | columnLower2[jColumn] = newLower; |
1595 | | columnUpper2[jColumn] = CoinMax(newLower, newUpper); |
1596 | | if (getRowStatus(iRow) != ClpSimplex::basic) { |
1597 | | if (getColumnStatus(iColumn) == ClpSimplex::basic) { |
1598 | | if (columnLower2[jColumn] == columnUpper2[jColumn]) { |
1599 | | // can only get here if will be fixed |
1600 | | small->setColumnStatus(jColumn, ClpSimplex::isFixed); |
1601 | | } else { |
1602 | | // solution is valid |
1603 | | if (fabs(columnActivity_[iColumn] - columnLower2[jColumn]) < |
1604 | | fabs(columnActivity_[iColumn] - columnUpper2[jColumn])) |
1605 | | small->setColumnStatus(jColumn, ClpSimplex::atLowerBound); |
1606 | | else |
1607 | | small->setColumnStatus(jColumn, ClpSimplex::atUpperBound); |
1608 | | } |
1609 | | } else { |
1610 | | //printf("what now neither basic\n"); |
1611 | | } |
1612 | | } |
1613 | | } |
1614 | | if (returnCode) { |
1615 | | delete small; |
1616 | | small = NULL; |
1617 | | } else if (tightenBounds && integerInformation) { |
1618 | | // See if we can tighten any bounds |
1619 | | // use rhs for upper and small duals for lower |
1620 | | double * up = rhs; |
1621 | | double * lo = small->dualRowSolution(); |
1622 | | const double * element = small->clpMatrix()->getElements(); |
1623 | | const int * row = small->clpMatrix()->getIndices(); |
1624 | | const CoinBigIndex * columnStart = small->clpMatrix()->getVectorStarts(); |
1625 | | //const int * columnLength = small->clpMatrix()->getVectorLengths(); |
1626 | | CoinZeroN(lo, numberRows2); |
1627 | | CoinZeroN(up, numberRows2); |
1628 | | for (int iColumn = 0; iColumn < numberColumns2; iColumn++) { |
1629 | | double upper = columnUpper2[iColumn]; |
1630 | | double lower = columnLower2[iColumn]; |
1631 | | //assert (columnLength[iColumn]==columnStart[iColumn+1]-columnStart[iColumn]); |
1632 | | for (CoinBigIndex j = columnStart[iColumn]; j < columnStart[iColumn+1]; j++) { |
1633 | | int iRow = row[j]; |
1634 | | double value = element[j]; |
1635 | | if (value > 0.0) { |
1636 | | if (upper < 1.0e20) |
1637 | | up[iRow] += upper * value; |
1638 | | else |
1639 | | up[iRow] = COIN_DBL_MAX; |
1640 | | if (lower > -1.0e20) |
1641 | | lo[iRow] += lower * value; |
1642 | | else |
1643 | | lo[iRow] = -COIN_DBL_MAX; |
1644 | | } else { |
1645 | | if (upper < 1.0e20) |
1646 | | lo[iRow] += upper * value; |
1647 | | else |
1648 | | lo[iRow] = -COIN_DBL_MAX; |
1649 | | if (lower > -1.0e20) |
1650 | | up[iRow] += lower * value; |
1651 | | else |
1652 | | up[iRow] = COIN_DBL_MAX; |
1653 | | } |
1654 | | } |
1655 | | } |
1656 | | double * rowLower2 = small->rowLower(); |
1657 | | double * rowUpper2 = small->rowUpper(); |
1658 | | bool feasible = true; |
1659 | | // make safer |
1660 | | for (int iRow = 0; iRow < numberRows2; iRow++) { |
1661 | | double lower = lo[iRow]; |
1662 | | if (lower > rowUpper2[iRow] + tolerance) { |
1663 | | feasible = false; |
1664 | | break; |
1665 | | } else { |
1666 | | lo[iRow] = CoinMin(lower - rowUpper2[iRow], 0.0) - tolerance; |
1667 | | } |
1668 | | double upper = up[iRow]; |
1669 | | if (upper < rowLower2[iRow] - tolerance) { |
1670 | | feasible = false; |
1671 | | break; |
1672 | | } else { |
1673 | | up[iRow] = CoinMax(upper - rowLower2[iRow], 0.0) + tolerance; |
1674 | | } |
1675 | | // tighten row bounds |
1676 | | if (lower>-1.0e10) |
1677 | | rowLower2[iRow] = CoinMax(rowLower2[iRow], |
1678 | | lower - 1.0e-6*(1.0+fabs(lower))); |
1679 | | if (upper<1.0e10) |
1680 | | rowUpper2[iRow] = CoinMin(rowUpper2[iRow], |
1681 | | upper + 1.0e-6*(1.0+fabs(upper))); |
1682 | | } |
1683 | | if (!feasible) { |
| 1736 | } else { |
| 1737 | // infeasible |
| 1738 | //printf("%d bounds infeasible %g %g tightened to %g %g\n", |
| 1739 | // iColumn,columnLower2[iColumn],columnUpper2[iColumn], |
| 1740 | // newLower,newUpper); |
| 1741 | #if 1 |
1686 | | } else { |
1687 | | // and tighten |
1688 | | for (int iColumn = 0; iColumn < numberColumns2; iColumn++) { |
1689 | | if (integerInformation[whichColumn[iColumn]]) { |
1690 | | double upper = columnUpper2[iColumn]; |
1691 | | double lower = columnLower2[iColumn]; |
1692 | | double newUpper = upper; |
1693 | | double newLower = lower; |
1694 | | double difference = upper - lower; |
1695 | | if (lower > -1000.0 && upper < 1000.0) { |
1696 | | for (CoinBigIndex j = columnStart[iColumn]; j < columnStart[iColumn+1]; j++) { |
1697 | | int iRow = row[j]; |
1698 | | double value = element[j]; |
1699 | | if (value > 0.0) { |
1700 | | double upWithOut = up[iRow] - value * difference; |
1701 | | if (upWithOut < 0.0) { |
1702 | | newLower = CoinMax(newLower, lower - (upWithOut + tolerance) / value); |
1703 | | } |
1704 | | double lowWithOut = lo[iRow] + value * difference; |
1705 | | if (lowWithOut > 0.0) { |
1706 | | newUpper = CoinMin(newUpper, upper - (lowWithOut - tolerance) / value); |
1707 | | } |
1708 | | } else { |
1709 | | double upWithOut = up[iRow] + value * difference; |
1710 | | if (upWithOut < 0.0) { |
1711 | | newUpper = CoinMin(newUpper, upper - (upWithOut + tolerance) / value); |
1712 | | } |
1713 | | double lowWithOut = lo[iRow] - value * difference; |
1714 | | if (lowWithOut > 0.0) { |
1715 | | newLower = CoinMax(newLower, lower - (lowWithOut - tolerance) / value); |
1716 | | } |
1717 | | } |
1718 | | } |
1719 | | if (newLower > lower || newUpper < upper) { |
1720 | | if (fabs(newUpper - floor(newUpper + 0.5)) > 1.0e-6) |
1721 | | newUpper = floor(newUpper); |
1722 | | else |
1723 | | newUpper = floor(newUpper + 0.5); |
1724 | | if (fabs(newLower - ceil(newLower - 0.5)) > 1.0e-6) |
1725 | | newLower = ceil(newLower); |
1726 | | else |
1727 | | newLower = ceil(newLower - 0.5); |
1728 | | // change may be too small - check |
1729 | | if (newLower > lower || newUpper < upper) { |
1730 | | if (newUpper >= newLower) { |
1731 | | // Could also tighten in this |
1732 | | //printf("%d bounds %g %g tightened to %g %g\n", |
1733 | | // iColumn,columnLower2[iColumn],columnUpper2[iColumn], |
1734 | | // newLower,newUpper); |
1735 | | #if 1 |
1736 | | columnUpper2[iColumn] = newUpper; |
1737 | | columnLower2[iColumn] = newLower; |
1738 | | columnUpper_[whichColumn[iColumn]] = newUpper; |
1739 | | columnLower_[whichColumn[iColumn]] = newLower; |
1740 | | #endif |
1741 | | // and adjust bounds on rows |
1742 | | newUpper -= upper; |
1743 | | newLower -= lower; |
1744 | | for (CoinBigIndex j = columnStart[iColumn]; j < columnStart[iColumn+1]; j++) { |
1745 | | int iRow = row[j]; |
1746 | | double value = element[j]; |
1747 | | if (value > 0.0) { |
1748 | | up[iRow] += newUpper * value; |
1749 | | lo[iRow] += newLower * value; |
1750 | | } else { |
1751 | | lo[iRow] += newUpper * value; |
1752 | | up[iRow] += newLower * value; |
1753 | | } |
1754 | | } |
1755 | | } else { |
1756 | | // infeasible |
1757 | | //printf("%d bounds infeasible %g %g tightened to %g %g\n", |
1758 | | // iColumn,columnLower2[iColumn],columnUpper2[iColumn], |
1759 | | // newLower,newUpper); |
1760 | | #if 1 |
1761 | | delete small; |
1762 | | small = NULL; |
1763 | | break; |
1764 | | #endif |
1765 | | } |
1766 | | } |
1767 | | } |
1768 | | } |
1769 | | } |
1770 | | } |
1771 | | } |
1772 | | } |
1773 | | } |
| 1744 | break; |
| 1745 | #endif |
| 1746 | } |
| 1747 | } |
| 1748 | } |
| 1749 | } |
| 1750 | } |
| 1751 | } |
| 1752 | } |
| 1753 | } |
| 1754 | } |
1817 | | for (int i = 0; i < small.numberRows(); i++) |
1818 | | assert (whichRow[i] >= 0 && whichRow[i] < numberRows_); |
1819 | | for (int i = 0; i < small.numberColumns(); i++) |
1820 | | assert (whichColumn[i] >= 0 && whichColumn[i] < numberColumns_); |
1821 | | #endif |
1822 | | getbackSolution(small, whichRow, whichColumn); |
1823 | | // and deal with status for bounds |
1824 | | const double * element = matrix_->getElements(); |
1825 | | const int * row = matrix_->getIndices(); |
1826 | | const CoinBigIndex * columnStart = matrix_->getVectorStarts(); |
1827 | | const int * columnLength = matrix_->getVectorLengths(); |
1828 | | double tolerance = primalTolerance(); |
1829 | | double djTolerance = dualTolerance(); |
1830 | | for (int jRow = nBound; jRow < 2 * numberRows_; jRow++) { |
1831 | | int iRow = whichRow[jRow]; |
1832 | | int iColumn = whichRow[jRow+numberRows_]; |
1833 | | if (getColumnStatus(iColumn) != ClpSimplex::basic) { |
1834 | | double lower = columnLower_[iColumn]; |
1835 | | double upper = columnUpper_[iColumn]; |
1836 | | double value = columnActivity_[iColumn]; |
1837 | | double djValue = reducedCost_[iColumn]; |
1838 | | dual_[iRow] = 0.0; |
1839 | | if (upper > lower) { |
1840 | | if (value < lower + tolerance && djValue > -djTolerance) { |
1841 | | setColumnStatus(iColumn, ClpSimplex::atLowerBound); |
1842 | | setRowStatus(iRow, ClpSimplex::basic); |
1843 | | } else if (value > upper - tolerance && djValue < djTolerance) { |
1844 | | setColumnStatus(iColumn, ClpSimplex::atUpperBound); |
1845 | | setRowStatus(iRow, ClpSimplex::basic); |
1846 | | } else { |
1847 | | // has to be basic |
1848 | | setColumnStatus(iColumn, ClpSimplex::basic); |
1849 | | reducedCost_[iColumn] = 0.0; |
1850 | | double value = 0.0; |
1851 | | for (CoinBigIndex j = columnStart[iColumn]; |
1852 | | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
1853 | | if (iRow == row[j]) { |
1854 | | value = element[j]; |
1855 | | break; |
1856 | | } |
1857 | | } |
1858 | | dual_[iRow] = djValue / value; |
1859 | | if (rowUpper_[iRow] > rowLower_[iRow]) { |
1860 | | if (fabs(rowActivity_[iRow] - rowLower_[iRow]) < |
1861 | | fabs(rowActivity_[iRow] - rowUpper_[iRow])) |
1862 | | setRowStatus(iRow, ClpSimplex::atLowerBound); |
1863 | | else |
1864 | | setRowStatus(iRow, ClpSimplex::atUpperBound); |
1865 | | } else { |
1866 | | setRowStatus(iRow, ClpSimplex::isFixed); |
1867 | | } |
1868 | | } |
1869 | | } else { |
1870 | | // row can always be basic |
1871 | | setRowStatus(iRow, ClpSimplex::basic); |
1872 | | } |
| 1797 | for (int i = 0; i < small.numberRows(); i++) |
| 1798 | assert(whichRow[i] >= 0 && whichRow[i] < numberRows_); |
| 1799 | for (int i = 0; i < small.numberColumns(); i++) |
| 1800 | assert(whichColumn[i] >= 0 && whichColumn[i] < numberColumns_); |
| 1801 | #endif |
| 1802 | getbackSolution(small, whichRow, whichColumn); |
| 1803 | // and deal with status for bounds |
| 1804 | const double *element = matrix_->getElements(); |
| 1805 | const int *row = matrix_->getIndices(); |
| 1806 | const CoinBigIndex *columnStart = matrix_->getVectorStarts(); |
| 1807 | const int *columnLength = matrix_->getVectorLengths(); |
| 1808 | double tolerance = primalTolerance(); |
| 1809 | double djTolerance = dualTolerance(); |
| 1810 | for (int jRow = nBound; jRow < 2 * numberRows_; jRow++) { |
| 1811 | int iRow = whichRow[jRow]; |
| 1812 | int iColumn = whichRow[jRow + numberRows_]; |
| 1813 | if (getColumnStatus(iColumn) != ClpSimplex::basic) { |
| 1814 | double lower = columnLower_[iColumn]; |
| 1815 | double upper = columnUpper_[iColumn]; |
| 1816 | double value = columnActivity_[iColumn]; |
| 1817 | double djValue = reducedCost_[iColumn]; |
| 1818 | dual_[iRow] = 0.0; |
| 1819 | if (upper > lower) { |
| 1820 | if (value < lower + tolerance && djValue > -djTolerance) { |
| 1821 | setColumnStatus(iColumn, ClpSimplex::atLowerBound); |
| 1822 | setRowStatus(iRow, ClpSimplex::basic); |
| 1823 | } else if (value > upper - tolerance && djValue < djTolerance) { |
| 1824 | setColumnStatus(iColumn, ClpSimplex::atUpperBound); |
| 1825 | setRowStatus(iRow, ClpSimplex::basic); |
| 1826 | } else { |
| 1827 | // has to be basic |
| 1828 | setColumnStatus(iColumn, ClpSimplex::basic); |
| 1829 | reducedCost_[iColumn] = 0.0; |
| 1830 | double value = 0.0; |
| 1831 | for (CoinBigIndex j = columnStart[iColumn]; |
| 1832 | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
| 1833 | if (iRow == row[j]) { |
| 1834 | value = element[j]; |
| 1835 | break; |
| 1836 | } |
| 1837 | } |
| 1838 | dual_[iRow] = djValue / value; |
| 1839 | if (rowUpper_[iRow] > rowLower_[iRow]) { |
| 1840 | if (fabs(rowActivity_[iRow] - rowLower_[iRow]) < fabs(rowActivity_[iRow] - rowUpper_[iRow])) |
| 1841 | setRowStatus(iRow, ClpSimplex::atLowerBound); |
| 1842 | else |
| 1843 | setRowStatus(iRow, ClpSimplex::atUpperBound); |
1898 | | // See if we can tighten any bounds |
1899 | | // use rhs for upper and small duals for lower |
1900 | | double * up = rhsSpace; |
1901 | | double * lo = dual_; |
1902 | | const double * element = matrix_->getElements(); |
1903 | | const int * row = matrix_->getIndices(); |
1904 | | const CoinBigIndex * columnStart = matrix_->getVectorStarts(); |
1905 | | const int * columnLength = matrix_->getVectorLengths(); |
1906 | | CoinZeroN(lo, numberRows_); |
1907 | | CoinZeroN(up, numberRows_); |
1908 | | for (int iColumn = 0; iColumn < numberColumns_; iColumn++) { |
1909 | | double upper = columnUpper_[iColumn]; |
1910 | | double lower = columnLower_[iColumn]; |
1911 | | //assert (columnLength[iColumn]==columnStart[iColumn+1]-columnStart[iColumn]); |
| 1876 | // See if we can tighten any bounds |
| 1877 | // use rhs for upper and small duals for lower |
| 1878 | double *up = rhsSpace; |
| 1879 | double *lo = dual_; |
| 1880 | const double *element = matrix_->getElements(); |
| 1881 | const int *row = matrix_->getIndices(); |
| 1882 | const CoinBigIndex *columnStart = matrix_->getVectorStarts(); |
| 1883 | const int *columnLength = matrix_->getVectorLengths(); |
| 1884 | CoinZeroN(lo, numberRows_); |
| 1885 | CoinZeroN(up, numberRows_); |
| 1886 | for (int iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 1887 | double upper = columnUpper_[iColumn]; |
| 1888 | double lower = columnLower_[iColumn]; |
| 1889 | //assert (columnLength[iColumn]==columnStart[iColumn+1]-columnStart[iColumn]); |
| 1890 | for (CoinBigIndex j = columnStart[iColumn]; |
| 1891 | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
| 1892 | int iRow = row[j]; |
| 1893 | double value = element[j]; |
| 1894 | if (value > 0.0) { |
| 1895 | if (upper < 1.0e20) |
| 1896 | up[iRow] += upper * value; |
| 1897 | else |
| 1898 | up[iRow] = COIN_DBL_MAX; |
| 1899 | if (lower > -1.0e20) |
| 1900 | lo[iRow] += lower * value; |
| 1901 | else |
| 1902 | lo[iRow] = -COIN_DBL_MAX; |
| 1903 | } else { |
| 1904 | if (upper < 1.0e20) |
| 1905 | lo[iRow] += upper * value; |
| 1906 | else |
| 1907 | lo[iRow] = -COIN_DBL_MAX; |
| 1908 | if (lower > -1.0e20) |
| 1909 | up[iRow] += lower * value; |
| 1910 | else |
| 1911 | up[iRow] = COIN_DBL_MAX; |
| 1912 | } |
| 1913 | } |
| 1914 | } |
| 1915 | bool feasible = true; |
| 1916 | // make safer |
| 1917 | double tolerance = primalTolerance(); |
| 1918 | for (int iRow = 0; iRow < numberRows_; iRow++) { |
| 1919 | double lower = lo[iRow]; |
| 1920 | if (lower > rowUpper_[iRow] + tolerance) { |
| 1921 | feasible = false; |
| 1922 | break; |
| 1923 | } else { |
| 1924 | lo[iRow] = CoinMin(lower - rowUpper_[iRow], 0.0) - tolerance; |
| 1925 | } |
| 1926 | double upper = up[iRow]; |
| 1927 | if (upper < rowLower_[iRow] - tolerance) { |
| 1928 | feasible = false; |
| 1929 | break; |
| 1930 | } else { |
| 1931 | up[iRow] = CoinMax(upper - rowLower_[iRow], 0.0) + tolerance; |
| 1932 | } |
| 1933 | } |
| 1934 | int numberTightened = 0; |
| 1935 | if (!feasible) { |
| 1936 | return -1; |
| 1937 | } else if (integerType_) { |
| 1938 | // and tighten |
| 1939 | for (int iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 1940 | if (integerType_[iColumn]) { |
| 1941 | double upper = columnUpper_[iColumn]; |
| 1942 | double lower = columnLower_[iColumn]; |
| 1943 | double newUpper = upper; |
| 1944 | double newLower = lower; |
| 1945 | double difference = upper - lower; |
| 1946 | if (lower > -1000.0 && upper < 1000.0) { |
1913 | | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
1914 | | int iRow = row[j]; |
1915 | | double value = element[j]; |
1916 | | if (value > 0.0) { |
1917 | | if (upper < 1.0e20) |
1918 | | up[iRow] += upper * value; |
1919 | | else |
1920 | | up[iRow] = COIN_DBL_MAX; |
1921 | | if (lower > -1.0e20) |
1922 | | lo[iRow] += lower * value; |
1923 | | else |
1924 | | lo[iRow] = -COIN_DBL_MAX; |
1925 | | } else { |
1926 | | if (upper < 1.0e20) |
1927 | | lo[iRow] += upper * value; |
1928 | | else |
1929 | | lo[iRow] = -COIN_DBL_MAX; |
1930 | | if (lower > -1.0e20) |
1931 | | up[iRow] += lower * value; |
1932 | | else |
1933 | | up[iRow] = COIN_DBL_MAX; |
1934 | | } |
1935 | | } |
1936 | | } |
1937 | | bool feasible = true; |
1938 | | // make safer |
1939 | | double tolerance = primalTolerance(); |
1940 | | for (int iRow = 0; iRow < numberRows_; iRow++) { |
1941 | | double lower = lo[iRow]; |
1942 | | if (lower > rowUpper_[iRow] + tolerance) { |
1943 | | feasible = false; |
1944 | | break; |
1945 | | } else { |
1946 | | lo[iRow] = CoinMin(lower - rowUpper_[iRow], 0.0) - tolerance; |
1947 | | } |
1948 | | double upper = up[iRow]; |
1949 | | if (upper < rowLower_[iRow] - tolerance) { |
1950 | | feasible = false; |
1951 | | break; |
1952 | | } else { |
1953 | | up[iRow] = CoinMax(upper - rowLower_[iRow], 0.0) + tolerance; |
1954 | | } |
1955 | | } |
1956 | | int numberTightened = 0; |
1957 | | if (!feasible) { |
1958 | | return -1; |
1959 | | } else if (integerType_) { |
1960 | | // and tighten |
1961 | | for (int iColumn = 0; iColumn < numberColumns_; iColumn++) { |
1962 | | if (integerType_[iColumn]) { |
1963 | | double upper = columnUpper_[iColumn]; |
1964 | | double lower = columnLower_[iColumn]; |
1965 | | double newUpper = upper; |
1966 | | double newLower = lower; |
1967 | | double difference = upper - lower; |
1968 | | if (lower > -1000.0 && upper < 1000.0) { |
1969 | | for (CoinBigIndex j = columnStart[iColumn]; |
1970 | | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
1971 | | int iRow = row[j]; |
1972 | | double value = element[j]; |
1973 | | if (value > 0.0) { |
1974 | | double upWithOut = up[iRow] - value * difference; |
1975 | | if (upWithOut < 0.0) { |
1976 | | newLower = CoinMax(newLower, lower - (upWithOut + tolerance) / value); |
1977 | | } |
1978 | | double lowWithOut = lo[iRow] + value * difference; |
1979 | | if (lowWithOut > 0.0) { |
1980 | | newUpper = CoinMin(newUpper, upper - (lowWithOut - tolerance) / value); |
1981 | | } |
1982 | | } else { |
1983 | | double upWithOut = up[iRow] + value * difference; |
1984 | | if (upWithOut < 0.0) { |
1985 | | newUpper = CoinMin(newUpper, upper - (upWithOut + tolerance) / value); |
1986 | | } |
1987 | | double lowWithOut = lo[iRow] - value * difference; |
1988 | | if (lowWithOut > 0.0) { |
1989 | | newLower = CoinMax(newLower, lower - (lowWithOut - tolerance) / value); |
1990 | | } |
1991 | | } |
1992 | | } |
1993 | | if (newLower > lower || newUpper < upper) { |
1994 | | if (fabs(newUpper - floor(newUpper + 0.5)) > 1.0e-6) |
1995 | | newUpper = floor(newUpper); |
1996 | | else |
1997 | | newUpper = floor(newUpper + 0.5); |
1998 | | if (fabs(newLower - ceil(newLower - 0.5)) > 1.0e-6) |
1999 | | newLower = ceil(newLower); |
2000 | | else |
2001 | | newLower = ceil(newLower - 0.5); |
2002 | | // change may be too small - check |
2003 | | if (newLower > lower || newUpper < upper) { |
2004 | | if (newUpper >= newLower) { |
2005 | | numberTightened++; |
2006 | | //printf("%d bounds %g %g tightened to %g %g\n", |
2007 | | // iColumn,columnLower_[iColumn],columnUpper_[iColumn], |
2008 | | // newLower,newUpper); |
2009 | | columnUpper_[iColumn] = newUpper; |
2010 | | columnLower_[iColumn] = newLower; |
2011 | | // and adjust bounds on rows |
2012 | | newUpper -= upper; |
2013 | | newLower -= lower; |
2014 | | for (CoinBigIndex j = columnStart[iColumn]; |
2015 | | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
2016 | | int iRow = row[j]; |
2017 | | double value = element[j]; |
2018 | | if (value > 0.0) { |
2019 | | up[iRow] += newUpper * value; |
2020 | | lo[iRow] += newLower * value; |
2021 | | } else { |
2022 | | lo[iRow] += newUpper * value; |
2023 | | up[iRow] += newLower * value; |
2024 | | } |
2025 | | } |
2026 | | } else { |
2027 | | // infeasible |
2028 | | //printf("%d bounds infeasible %g %g tightened to %g %g\n", |
2029 | | // iColumn,columnLower_[iColumn],columnUpper_[iColumn], |
2030 | | // newLower,newUpper); |
2031 | | return -1; |
2032 | | } |
2033 | | } |
2034 | | } |
2035 | | } |
2036 | | } |
2037 | | } |
2038 | | } |
2039 | | return numberTightened; |
| 1948 | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
| 1949 | int iRow = row[j]; |
| 1950 | double value = element[j]; |
| 1951 | if (value > 0.0) { |
| 1952 | double upWithOut = up[iRow] - value * difference; |
| 1953 | if (upWithOut < 0.0) { |
| 1954 | newLower = CoinMax(newLower, lower - (upWithOut + tolerance) / value); |
| 1955 | } |
| 1956 | double lowWithOut = lo[iRow] + value * difference; |
| 1957 | if (lowWithOut > 0.0) { |
| 1958 | newUpper = CoinMin(newUpper, upper - (lowWithOut - tolerance) / value); |
| 1959 | } |
| 1960 | } else { |
| 1961 | double upWithOut = up[iRow] + value * difference; |
| 1962 | if (upWithOut < 0.0) { |
| 1963 | newUpper = CoinMin(newUpper, upper - (upWithOut + tolerance) / value); |
| 1964 | } |
| 1965 | double lowWithOut = lo[iRow] - value * difference; |
| 1966 | if (lowWithOut > 0.0) { |
| 1967 | newLower = CoinMax(newLower, lower - (lowWithOut - tolerance) / value); |
| 1968 | } |
| 1969 | } |
| 1970 | } |
| 1971 | if (newLower > lower || newUpper < upper) { |
| 1972 | if (fabs(newUpper - floor(newUpper + 0.5)) > 1.0e-6) |
| 1973 | newUpper = floor(newUpper); |
| 1974 | else |
| 1975 | newUpper = floor(newUpper + 0.5); |
| 1976 | if (fabs(newLower - ceil(newLower - 0.5)) > 1.0e-6) |
| 1977 | newLower = ceil(newLower); |
| 1978 | else |
| 1979 | newLower = ceil(newLower - 0.5); |
| 1980 | // change may be too small - check |
| 1981 | if (newLower > lower || newUpper < upper) { |
| 1982 | if (newUpper >= newLower) { |
| 1983 | numberTightened++; |
| 1984 | //printf("%d bounds %g %g tightened to %g %g\n", |
| 1985 | // iColumn,columnLower_[iColumn],columnUpper_[iColumn], |
| 1986 | // newLower,newUpper); |
| 1987 | columnUpper_[iColumn] = newUpper; |
| 1988 | columnLower_[iColumn] = newLower; |
| 1989 | // and adjust bounds on rows |
| 1990 | newUpper -= upper; |
| 1991 | newLower -= lower; |
| 1992 | for (CoinBigIndex j = columnStart[iColumn]; |
| 1993 | j < columnStart[iColumn] + columnLength[iColumn]; j++) { |
| 1994 | int iRow = row[j]; |
| 1995 | double value = element[j]; |
| 1996 | if (value > 0.0) { |
| 1997 | up[iRow] += newUpper * value; |
| 1998 | lo[iRow] += newLower * value; |
| 1999 | } else { |
| 2000 | lo[iRow] += newUpper * value; |
| 2001 | up[iRow] += newLower * value; |
| 2002 | } |
| 2003 | } |
| 2004 | } else { |
| 2005 | // infeasible |
| 2006 | //printf("%d bounds infeasible %g %g tightened to %g %g\n", |
| 2007 | // iColumn,columnLower_[iColumn],columnUpper_[iColumn], |
| 2008 | // newLower,newUpper); |
| 2009 | return -1; |
| 2010 | } |
| 2011 | } |
| 2012 | } |
| 2013 | } |
| 2014 | } |
| 2015 | } |
| 2016 | } |
| 2017 | return numberTightened; |
| 2059 | if (!returnCode) { |
| 2060 | // Find theta when bounds will cross over and create arrays |
| 2061 | int numberTotal = numberRows_ + numberColumns_; |
| 2062 | chgLower = new double[numberTotal]; |
| 2063 | memset(chgLower, 0, numberTotal * sizeof(double)); |
| 2064 | chgUpper = new double[numberTotal]; |
| 2065 | memset(chgUpper, 0, numberTotal * sizeof(double)); |
| 2066 | chgObjective = new double[numberTotal]; |
| 2067 | memset(chgObjective, 0, numberTotal * sizeof(double)); |
| 2068 | assert(!rowScale_); |
| 2069 | double maxTheta = 1.0e50; |
| 2070 | if (lowerChangeRhs || upperChangeRhs) { |
| 2071 | for (iRow = 0; iRow < numberRows_; iRow++) { |
| 2072 | double lower = rowLower_[iRow]; |
| 2073 | double upper = rowUpper_[iRow]; |
| 2074 | if (lower > upper) { |
| 2075 | maxTheta = -1.0; |
| 2076 | break; |
| 2077 | } |
| 2078 | double lowerChange = (lowerChangeRhs) ? lowerChangeRhs[iRow] : 0.0; |
| 2079 | double upperChange = (upperChangeRhs) ? upperChangeRhs[iRow] : 0.0; |
| 2080 | if (lower > -1.0e20 && upper < 1.0e20) { |
| 2081 | if (lower + maxTheta * lowerChange > upper + maxTheta * upperChange) { |
| 2082 | maxTheta = (upper - lower) / (lowerChange - upperChange); |
| 2083 | } |
| 2084 | } |
| 2085 | if (lower > -1.0e20) { |
| 2086 | lower_[numberColumns_ + iRow] += startingTheta * lowerChange; |
| 2087 | chgLower[numberColumns_ + iRow] = lowerChange; |
| 2088 | } |
| 2089 | if (upper < 1.0e20) { |
| 2090 | upper_[numberColumns_ + iRow] += startingTheta * upperChange; |
| 2091 | chgUpper[numberColumns_ + iRow] = upperChange; |
| 2092 | } |
| 2093 | } |
| 2094 | } |
| 2095 | if (maxTheta > 0.0) { |
| 2096 | if (lowerChangeBound || upperChangeBound) { |
| 2097 | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 2098 | double lower = columnLower_[iColumn]; |
| 2099 | double upper = columnUpper_[iColumn]; |
| 2100 | if (lower > upper) { |
| 2101 | maxTheta = -1.0; |
| 2102 | break; |
| 2103 | } |
| 2104 | double lowerChange = (lowerChangeBound) ? lowerChangeBound[iColumn] : 0.0; |
| 2105 | double upperChange = (upperChangeBound) ? upperChangeBound[iColumn] : 0.0; |
| 2106 | if (lower > -1.0e20 && upper < 1.0e20) { |
| 2107 | if (lower + maxTheta * lowerChange > upper + maxTheta * upperChange) { |
| 2108 | maxTheta = (upper - lower) / (lowerChange - upperChange); |
| 2109 | } |
| 2110 | } |
| 2111 | if (lower > -1.0e20) { |
| 2112 | lower_[iColumn] += startingTheta * lowerChange; |
| 2113 | chgLower[iColumn] = lowerChange; |
| 2114 | } |
| 2115 | if (upper < 1.0e20) { |
| 2116 | upper_[iColumn] += startingTheta * upperChange; |
| 2117 | chgUpper[iColumn] = upperChange; |
| 2118 | } |
| 2119 | } |
| 2120 | } |
| 2121 | if (maxTheta == 1.0e50) |
| 2122 | maxTheta = COIN_DBL_MAX; |
| 2123 | } |
| 2124 | if (maxTheta < 0.0) { |
| 2125 | // bad ranges or initial |
| 2126 | returnCode = -1; |
| 2127 | } |
| 2128 | if (maxTheta < endingTheta) { |
| 2129 | char line[100]; |
| 2130 | sprintf(line, "Crossover considerations reduce ending theta from %g to %g\n", |
| 2131 | endingTheta, maxTheta); |
| 2132 | handler_->message(CLP_GENERAL, messages_) |
| 2133 | << line << CoinMessageEol; |
| 2134 | endingTheta = maxTheta; |
| 2135 | } |
| 2136 | if (endingTheta < startingTheta) { |
| 2137 | // bad initial |
| 2138 | returnCode = -2; |
| 2139 | } |
| 2140 | } |
| 2141 | double saveEndingTheta = endingTheta; |
| 2142 | if (!returnCode) { |
| 2143 | if (changeObjective) { |
| 2144 | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 2145 | chgObjective[iColumn] = changeObjective[iColumn]; |
| 2146 | cost_[iColumn] += startingTheta * changeObjective[iColumn]; |
| 2147 | } |
| 2148 | } |
| 2149 | double *saveDuals = NULL; |
| 2150 | reinterpret_cast< ClpSimplexDual * >(this)->gutsOfDual(0, saveDuals, -1, data); |
| 2151 | assert(!problemStatus_); |
| 2152 | for (int i = 0; i < numberRows_ + numberColumns_; i++) |
| 2153 | setFakeBound(i, noFake); |
| 2154 | // Now do parametrics |
| 2155 | handler_->message(CLP_PARAMETRICS_STATS, messages_) |
| 2156 | << startingTheta << objectiveValue() << CoinMessageEol; |
| 2157 | while (!returnCode) { |
| 2158 | //assert (reportIncrement); |
| 2159 | parametricsData paramData; |
| 2160 | paramData.startingTheta = startingTheta; |
| 2161 | paramData.endingTheta = endingTheta; |
| 2162 | paramData.maxTheta = COIN_DBL_MAX; |
| 2163 | paramData.lowerChange = chgLower; |
| 2164 | paramData.upperChange = chgUpper; |
| 2165 | returnCode = parametricsLoop(paramData, reportIncrement, |
| 2166 | chgLower, chgUpper, chgObjective, data, |
| 2167 | canTryQuick); |
| 2168 | startingTheta = paramData.startingTheta; |
| 2169 | endingTheta = paramData.endingTheta; |
| 2170 | if (!returnCode) { |
| 2171 | //double change = endingTheta-startingTheta; |
| 2172 | startingTheta = endingTheta; |
| 2173 | endingTheta = saveEndingTheta; |
| 2174 | //for (int i=0;i<numberTotal;i++) { |
| 2175 | //lower_[i] += change*chgLower[i]; |
| 2176 | //upper_[i] += change*chgUpper[i]; |
| 2177 | //cost_[i] += change*chgObjective[i]; |
| 2178 | //} |
| 2179 | handler_->message(CLP_PARAMETRICS_STATS, messages_) |
| 2180 | << startingTheta << objectiveValue() << CoinMessageEol; |
| 2181 | if (startingTheta >= endingTheta) |
| 2182 | break; |
| 2183 | } else if (returnCode == -1) { |
| 2184 | // trouble - do external solve |
| 2185 | needToDoSomething = true; |
| 2186 | } else if (problemStatus_ == 1) { |
| 2187 | // can't move any further |
| 2188 | if (!canTryQuick) { |
| 2189 | handler_->message(CLP_PARAMETRICS_STATS, messages_) |
| 2190 | << endingTheta << objectiveValue() << CoinMessageEol; |
| 2191 | problemStatus_ = 0; |
| 2192 | } |
| 2193 | } else { |
| 2194 | abort(); |
| 2195 | } |
| 2196 | } |
| 2197 | } |
| 2198 | reinterpret_cast< ClpSimplexDual * >(this)->finishSolve(0); |
2083 | | if (!returnCode) { |
2084 | | // Find theta when bounds will cross over and create arrays |
2085 | | int numberTotal = numberRows_ + numberColumns_; |
2086 | | chgLower = new double[numberTotal]; |
2087 | | memset(chgLower, 0, numberTotal * sizeof(double)); |
2088 | | chgUpper = new double[numberTotal]; |
2089 | | memset(chgUpper, 0, numberTotal * sizeof(double)); |
2090 | | chgObjective = new double[numberTotal]; |
2091 | | memset(chgObjective, 0, numberTotal * sizeof(double)); |
2092 | | assert (!rowScale_); |
2093 | | double maxTheta = 1.0e50; |
2094 | | if (lowerChangeRhs || upperChangeRhs) { |
2095 | | for (iRow = 0; iRow < numberRows_; iRow++) { |
2096 | | double lower = rowLower_[iRow]; |
2097 | | double upper = rowUpper_[iRow]; |
2098 | | if (lower > upper) { |
2099 | | maxTheta = -1.0; |
2100 | | break; |
2101 | | } |
2102 | | double lowerChange = (lowerChangeRhs) ? lowerChangeRhs[iRow] : 0.0; |
2103 | | double upperChange = (upperChangeRhs) ? upperChangeRhs[iRow] : 0.0; |
2104 | | if (lower > -1.0e20 && upper < 1.0e20) { |
2105 | | if (lower + maxTheta * lowerChange > upper + maxTheta * upperChange) { |
2106 | | maxTheta = (upper - lower) / (lowerChange - upperChange); |
2107 | | } |
2108 | | } |
2109 | | if (lower > -1.0e20) { |
2110 | | lower_[numberColumns_+iRow] += startingTheta * lowerChange; |
2111 | | chgLower[numberColumns_+iRow] = lowerChange; |
2112 | | } |
2113 | | if (upper < 1.0e20) { |
2114 | | upper_[numberColumns_+iRow] += startingTheta * upperChange; |
2115 | | chgUpper[numberColumns_+iRow] = upperChange; |
2116 | | } |
2117 | | } |
2118 | | } |
2119 | | if (maxTheta > 0.0) { |
2120 | | if (lowerChangeBound || upperChangeBound) { |
2121 | | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
2122 | | double lower = columnLower_[iColumn]; |
2123 | | double upper = columnUpper_[iColumn]; |
2124 | | if (lower > upper) { |
2125 | | maxTheta = -1.0; |
2126 | | break; |
2127 | | } |
2128 | | double lowerChange = (lowerChangeBound) ? lowerChangeBound[iColumn] : 0.0; |
2129 | | double upperChange = (upperChangeBound) ? upperChangeBound[iColumn] : 0.0; |
2130 | | if (lower > -1.0e20 && upper < 1.0e20) { |
2131 | | if (lower + maxTheta * lowerChange > upper + maxTheta * upperChange) { |
2132 | | maxTheta = (upper - lower) / (lowerChange - upperChange); |
2133 | | } |
2134 | | } |
2135 | | if (lower > -1.0e20) { |
2136 | | lower_[iColumn] += startingTheta * lowerChange; |
2137 | | chgLower[iColumn] = lowerChange; |
2138 | | } |
2139 | | if (upper < 1.0e20) { |
2140 | | upper_[iColumn] += startingTheta * upperChange; |
2141 | | chgUpper[iColumn] = upperChange; |
2142 | | } |
2143 | | } |
2144 | | } |
2145 | | if (maxTheta == 1.0e50) |
2146 | | maxTheta = COIN_DBL_MAX; |
2147 | | } |
2148 | | if (maxTheta < 0.0) { |
2149 | | // bad ranges or initial |
2150 | | returnCode = -1; |
2151 | | } |
2152 | | if (maxTheta < endingTheta) { |
2153 | | char line[100]; |
2154 | | sprintf(line,"Crossover considerations reduce ending theta from %g to %g\n", |
2155 | | endingTheta,maxTheta); |
2156 | | handler_->message(CLP_GENERAL,messages_) |
2157 | | << line << CoinMessageEol; |
2158 | | endingTheta = maxTheta; |
2159 | | } |
2160 | | if (endingTheta < startingTheta) { |
2161 | | // bad initial |
2162 | | returnCode = -2; |
2163 | | } |
2164 | | } |
2165 | | double saveEndingTheta = endingTheta; |
2166 | | if (!returnCode) { |
2167 | | if (changeObjective) { |
2168 | | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
2169 | | chgObjective[iColumn] = changeObjective[iColumn]; |
2170 | | cost_[iColumn] += startingTheta * changeObjective[iColumn]; |
2171 | | } |
2172 | | } |
2173 | | double * saveDuals = NULL; |
2174 | | reinterpret_cast<ClpSimplexDual *> (this)->gutsOfDual(0, saveDuals, -1, data); |
2175 | | assert (!problemStatus_); |
2176 | | for (int i=0;i<numberRows_+numberColumns_;i++) |
2177 | | setFakeBound(i, noFake); |
2178 | | // Now do parametrics |
2179 | | handler_->message(CLP_PARAMETRICS_STATS, messages_) |
2180 | | << startingTheta << objectiveValue() << CoinMessageEol; |
2181 | | while (!returnCode) { |
2182 | | //assert (reportIncrement); |
2183 | | parametricsData paramData; |
2184 | | paramData.startingTheta=startingTheta; |
2185 | | paramData.endingTheta=endingTheta; |
2186 | | paramData.maxTheta=COIN_DBL_MAX; |
2187 | | paramData.lowerChange = chgLower; |
2188 | | paramData.upperChange = chgUpper; |
2189 | | returnCode = parametricsLoop(paramData, reportIncrement, |
2190 | | chgLower, chgUpper, chgObjective, data, |
2191 | | canTryQuick); |
2192 | | startingTheta=paramData.startingTheta; |
2193 | | endingTheta=paramData.endingTheta; |
2194 | | if (!returnCode) { |
2195 | | //double change = endingTheta-startingTheta; |
2196 | | startingTheta = endingTheta; |
2197 | | endingTheta = saveEndingTheta; |
2198 | | //for (int i=0;i<numberTotal;i++) { |
2199 | | //lower_[i] += change*chgLower[i]; |
2200 | | //upper_[i] += change*chgUpper[i]; |
2201 | | //cost_[i] += change*chgObjective[i]; |
2202 | | //} |
2203 | | handler_->message(CLP_PARAMETRICS_STATS, messages_) |
2204 | | << startingTheta << objectiveValue() << CoinMessageEol; |
2205 | | if (startingTheta >= endingTheta) |
2206 | | break; |
2207 | | } else if (returnCode == -1) { |
2208 | | // trouble - do external solve |
2209 | | needToDoSomething = true; |
2210 | | } else if (problemStatus_==1) { |
2211 | | // can't move any further |
2212 | | if (!canTryQuick) { |
2213 | | handler_->message(CLP_PARAMETRICS_STATS, messages_) |
2214 | | << endingTheta << objectiveValue() << CoinMessageEol; |
2215 | | problemStatus_=0; |
2216 | | } |
2217 | | } else { |
2218 | | abort(); |
2219 | | } |
2220 | | } |
2221 | | } |
2222 | | reinterpret_cast<ClpSimplexDual *> (this)->finishSolve(0); |
2223 | | |
2224 | | delete dualRowPivot_; |
2225 | | dualRowPivot_ = savePivot; |
2226 | | // Restore any saved stuff |
2227 | | restoreData(data); |
2228 | | if (needToDoSomething) { |
2229 | | double saveStartingTheta = startingTheta; // known to be feasible |
2230 | | int cleanedUp = 1; |
2231 | | while (cleanedUp) { |
2232 | | // tweak |
2233 | | if (cleanedUp == 1) { |
2234 | | if (!reportIncrement) |
2235 | | startingTheta = CoinMin(startingTheta + 1.0e-5, saveEndingTheta); |
2236 | | else |
2237 | | startingTheta = CoinMin(startingTheta + reportIncrement, saveEndingTheta); |
2238 | | } else { |
2239 | | // restoring to go slowly |
2240 | | startingTheta = saveStartingTheta; |
2241 | | } |
2242 | | // only works if not scaled |
2243 | | int i; |
2244 | | const double * obj1 = objective(); |
2245 | | double * obj2 = copyModel.objective(); |
2246 | | const double * lower1 = columnLower_; |
2247 | | double * lower2 = copyModel.columnLower(); |
2248 | | const double * upper1 = columnUpper_; |
2249 | | double * upper2 = copyModel.columnUpper(); |
2250 | | for (i = 0; i < numberColumns_; i++) { |
2251 | | obj2[i] = obj1[i] + startingTheta * chgObjective[i]; |
2252 | | lower2[i] = lower1[i] + startingTheta * chgLower[i]; |
2253 | | upper2[i] = upper1[i] + startingTheta * chgUpper[i]; |
2254 | | } |
2255 | | lower1 = rowLower_; |
2256 | | lower2 = copyModel.rowLower(); |
2257 | | upper1 = rowUpper_; |
2258 | | upper2 = copyModel.rowUpper(); |
2259 | | for (i = 0; i < numberRows_; i++) { |
2260 | | lower2[i] = lower1[i] + startingTheta * chgLower[i+numberColumns_]; |
2261 | | upper2[i] = upper1[i] + startingTheta * chgUpper[i+numberColumns_]; |
2262 | | } |
2263 | | copyModel.dual(); |
2264 | | if (copyModel.problemStatus()) { |
2265 | | char line[100]; |
2266 | | sprintf(line,"Can not get to theta of %g\n", startingTheta); |
2267 | | handler_->message(CLP_GENERAL,messages_) |
2268 | | << line << CoinMessageEol; |
2269 | | canTryQuick = false; // do slowly to get exact amount |
2270 | | // back to last known good |
2271 | | if (cleanedUp == 1) |
2272 | | cleanedUp = 2; |
2273 | | else |
2274 | | abort(); |
2275 | | } else { |
2276 | | // and move stuff back |
2277 | | int numberTotal = numberRows_ + numberColumns_; |
2278 | | CoinMemcpyN(copyModel.statusArray(), numberTotal, status_); |
2279 | | CoinMemcpyN(copyModel.primalColumnSolution(), numberColumns_, columnActivity_); |
2280 | | CoinMemcpyN(copyModel.primalRowSolution(), numberRows_, rowActivity_); |
2281 | | cleanedUp = 0; |
2282 | | } |
2283 | | } |
2284 | | } |
2285 | | delete [] chgLower; |
2286 | | delete [] chgUpper; |
2287 | | delete [] chgObjective; |
2288 | | } |
2289 | | perturbation_ = savePerturbation; |
2290 | | char line[100]; |
2291 | | sprintf(line,"Ending theta %g\n", endingTheta); |
2292 | | handler_->message(CLP_GENERAL,messages_) |
2293 | | << line << CoinMessageEol; |
2294 | | return problemStatus_; |
| 2200 | delete dualRowPivot_; |
| 2201 | dualRowPivot_ = savePivot; |
| 2202 | // Restore any saved stuff |
| 2203 | restoreData(data); |
| 2204 | if (needToDoSomething) { |
| 2205 | double saveStartingTheta = startingTheta; // known to be feasible |
| 2206 | int cleanedUp = 1; |
| 2207 | while (cleanedUp) { |
| 2208 | // tweak |
| 2209 | if (cleanedUp == 1) { |
| 2210 | if (!reportIncrement) |
| 2211 | startingTheta = CoinMin(startingTheta + 1.0e-5, saveEndingTheta); |
| 2212 | else |
| 2213 | startingTheta = CoinMin(startingTheta + reportIncrement, saveEndingTheta); |
| 2214 | } else { |
| 2215 | // restoring to go slowly |
| 2216 | startingTheta = saveStartingTheta; |
| 2217 | } |
| 2218 | // only works if not scaled |
| 2219 | int i; |
| 2220 | const double *obj1 = objective(); |
| 2221 | double *obj2 = copyModel.objective(); |
| 2222 | const double *lower1 = columnLower_; |
| 2223 | double *lower2 = copyModel.columnLower(); |
| 2224 | const double *upper1 = columnUpper_; |
| 2225 | double *upper2 = copyModel.columnUpper(); |
| 2226 | for (i = 0; i < numberColumns_; i++) { |
| 2227 | obj2[i] = obj1[i] + startingTheta * chgObjective[i]; |
| 2228 | lower2[i] = lower1[i] + startingTheta * chgLower[i]; |
| 2229 | upper2[i] = upper1[i] + startingTheta * chgUpper[i]; |
| 2230 | } |
| 2231 | lower1 = rowLower_; |
| 2232 | lower2 = copyModel.rowLower(); |
| 2233 | upper1 = rowUpper_; |
| 2234 | upper2 = copyModel.rowUpper(); |
| 2235 | for (i = 0; i < numberRows_; i++) { |
| 2236 | lower2[i] = lower1[i] + startingTheta * chgLower[i + numberColumns_]; |
| 2237 | upper2[i] = upper1[i] + startingTheta * chgUpper[i + numberColumns_]; |
| 2238 | } |
| 2239 | copyModel.dual(); |
| 2240 | if (copyModel.problemStatus()) { |
| 2241 | char line[100]; |
| 2242 | sprintf(line, "Can not get to theta of %g\n", startingTheta); |
| 2243 | handler_->message(CLP_GENERAL, messages_) |
| 2244 | << line << CoinMessageEol; |
| 2245 | canTryQuick = false; // do slowly to get exact amount |
| 2246 | // back to last known good |
| 2247 | if (cleanedUp == 1) |
| 2248 | cleanedUp = 2; |
| 2249 | else |
| 2250 | abort(); |
| 2251 | } else { |
| 2252 | // and move stuff back |
| 2253 | int numberTotal = numberRows_ + numberColumns_; |
| 2254 | CoinMemcpyN(copyModel.statusArray(), numberTotal, status_); |
| 2255 | CoinMemcpyN(copyModel.primalColumnSolution(), numberColumns_, columnActivity_); |
| 2256 | CoinMemcpyN(copyModel.primalRowSolution(), numberRows_, rowActivity_); |
| 2257 | cleanedUp = 0; |
| 2258 | } |
| 2259 | } |
| 2260 | } |
| 2261 | delete[] chgLower; |
| 2262 | delete[] chgUpper; |
| 2263 | delete[] chgObjective; |
| 2264 | } |
| 2265 | perturbation_ = savePerturbation; |
| 2266 | char line[100]; |
| 2267 | sprintf(line, "Ending theta %g\n", endingTheta); |
| 2268 | handler_->message(CLP_GENERAL, messages_) |
| 2269 | << line << CoinMessageEol; |
| 2270 | return problemStatus_; |
2483 | | if (!strncmp(line, "ENDATA", 6)|| |
2484 | | !strncmp(line, "COLUMN",6)) |
2485 | | break; |
2486 | | nLine++; |
2487 | | iRow = -1; |
2488 | | double upper = 0.0; |
2489 | | double lower = 0.0; |
2490 | | char * pos = line; |
2491 | | char * put = line; |
2492 | | while (*pos >= ' ' && *pos != '\n') { |
2493 | | if (*pos != ' ' && *pos != '\t') { |
2494 | | *put = *pos; |
2495 | | put++; |
2496 | | } |
2497 | | pos++; |
2498 | | } |
2499 | | *put = '\0'; |
2500 | | pos = line; |
2501 | | for (int i = 0; i < nAcross; i++) { |
2502 | | char * comma = strchr(pos, ','); |
2503 | | if (comma) { |
2504 | | *comma = '\0'; |
2505 | | } else if (i < nAcross - 1) { |
2506 | | nBadLine++; |
2507 | | break; |
2508 | | } |
2509 | | switch (orderRow[i]) { |
2510 | | // name |
2511 | | case 0: |
2512 | | // For large problems this could be slow |
2513 | | for (iRow = 0; iRow < numberRows_; iRow++) { |
2514 | | if (!strcmp(rowNames[iRow], pos)) |
2515 | | break; |
2516 | | } |
2517 | | if (iRow == numberRows_) |
2518 | | iRow = -1; |
2519 | | break; |
2520 | | // number |
2521 | | case 1: |
2522 | | iRow = atoi(pos); |
2523 | | if (iRow < 0 || iRow >= numberRows_) |
2524 | | iRow = -1; |
2525 | | break; |
2526 | | // lower |
2527 | | case 2: |
2528 | | upper = atof(pos); |
2529 | | break; |
2530 | | // upper |
2531 | | case 3: |
2532 | | lower = atof(pos); |
2533 | | break; |
2534 | | // rhs |
2535 | | case 4: |
2536 | | lower = atof(pos); |
2537 | | upper = lower; |
2538 | | break; |
2539 | | } |
2540 | | if (comma) { |
2541 | | *comma = ','; |
2542 | | pos = comma + 1; |
2543 | | } |
2544 | | } |
2545 | | if (iRow >= 0) { |
2546 | | if (rowLower_[iRow]>-1.0e20) |
2547 | | lowerRowMove[iRow] = lower; |
2548 | | else |
2549 | | lowerRowMove[iRow]=0.0; |
2550 | | if (rowUpper_[iRow]<1.0e20) |
2551 | | upperRowMove[iRow] = upper; |
2552 | | else |
2553 | | upperRowMove[iRow] = lower; |
2554 | | } else { |
2555 | | nBadName++; |
2556 | | if(saveLine[0]=='\0') |
2557 | | strcpy(saveLine,line); |
2558 | | } |
2559 | | } |
2560 | | sprintf(line,"%d Row fields and %d records", nAcross, nLine); |
2561 | | handler_->message(CLP_GENERAL,messages_) |
2562 | | << line << CoinMessageEol; |
| 2455 | if (!strncmp(line, "ENDATA", 6) || !strncmp(line, "COLUMN", 6)) |
| 2456 | break; |
| 2457 | nLine++; |
| 2458 | iRow = -1; |
| 2459 | double upper = 0.0; |
| 2460 | double lower = 0.0; |
| 2461 | char *pos = line; |
| 2462 | char *put = line; |
| 2463 | while (*pos >= ' ' && *pos != '\n') { |
| 2464 | if (*pos != ' ' && *pos != '\t') { |
| 2465 | *put = *pos; |
| 2466 | put++; |
| 2467 | } |
| 2468 | pos++; |
| 2469 | } |
| 2470 | *put = '\0'; |
| 2471 | pos = line; |
| 2472 | for (int i = 0; i < nAcross; i++) { |
| 2473 | char *comma = strchr(pos, ','); |
| 2474 | if (comma) { |
| 2475 | *comma = '\0'; |
| 2476 | } else if (i < nAcross - 1) { |
| 2477 | nBadLine++; |
| 2478 | break; |
| 2479 | } |
| 2480 | switch (orderRow[i]) { |
| 2481 | // name |
| 2482 | case 0: |
| 2483 | // For large problems this could be slow |
| 2484 | for (iRow = 0; iRow < numberRows_; iRow++) { |
| 2485 | if (!strcmp(rowNames[iRow], pos)) |
| 2486 | break; |
| 2487 | } |
| 2488 | if (iRow == numberRows_) |
| 2489 | iRow = -1; |
| 2490 | break; |
| 2491 | // number |
| 2492 | case 1: |
| 2493 | iRow = atoi(pos); |
| 2494 | if (iRow < 0 || iRow >= numberRows_) |
| 2495 | iRow = -1; |
| 2496 | break; |
| 2497 | // lower |
| 2498 | case 2: |
| 2499 | upper = atof(pos); |
| 2500 | break; |
| 2501 | // upper |
| 2502 | case 3: |
| 2503 | lower = atof(pos); |
| 2504 | break; |
| 2505 | // rhs |
| 2506 | case 4: |
| 2507 | lower = atof(pos); |
| 2508 | upper = lower; |
| 2509 | break; |
| 2510 | } |
| 2511 | if (comma) { |
| 2512 | *comma = ','; |
| 2513 | pos = comma + 1; |
| 2514 | } |
| 2515 | } |
| 2516 | if (iRow >= 0) { |
| 2517 | if (rowLower_[iRow] > -1.0e20) |
| 2518 | lowerRowMove[iRow] = lower; |
| 2519 | else |
| 2520 | lowerRowMove[iRow] = 0.0; |
| 2521 | if (rowUpper_[iRow] < 1.0e20) |
| 2522 | upperRowMove[iRow] = upper; |
| 2523 | else |
| 2524 | upperRowMove[iRow] = lower; |
| 2525 | } else { |
| 2526 | nBadName++; |
| 2527 | if (saveLine[0] == '\0') |
| 2528 | strcpy(saveLine, line); |
| 2529 | } |
| 2530 | } |
| 2531 | sprintf(line, "%d Row fields and %d records", nAcross, nLine); |
| 2532 | handler_->message(CLP_GENERAL, messages_) |
| 2533 | << line << CoinMessageEol; |
2641 | | char ** columnNames = new char * [numberColumns_]; |
2642 | | int iColumn; |
2643 | | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
2644 | | columnNames[iColumn] = |
2645 | | CoinStrdup(columnName(iColumn).c_str()); |
2646 | | } |
2647 | | lowerColumnMove = new double [numberColumns_]; |
2648 | | memset(lowerColumnMove,0,numberColumns_*sizeof(double)); |
2649 | | upperColumnMove = new double [numberColumns_]; |
2650 | | memset(upperColumnMove,0,numberColumns_*sizeof(double)); |
2651 | | objectiveMove = new double [numberColumns_]; |
2652 | | memset(objectiveMove,0,numberColumns_*sizeof(double)); |
2653 | | int nLine = 0; |
2654 | | int nBadLine = 0; |
2655 | | int nBadName = 0; |
2656 | | while (fgets(line, 200, fp)) { |
2657 | | if (!strncmp(line, "ENDATA", 6)) |
2658 | | break; |
2659 | | nLine++; |
2660 | | iColumn = -1; |
2661 | | double upper = 0.0; |
2662 | | double lower = 0.0; |
2663 | | double obj =0.0; |
2664 | | char * pos = line; |
2665 | | char * put = line; |
2666 | | while (*pos >= ' ' && *pos != '\n') { |
2667 | | if (*pos != ' ' && *pos != '\t') { |
2668 | | *put = *pos; |
2669 | | put++; |
2670 | | } |
2671 | | pos++; |
2672 | | } |
2673 | | *put = '\0'; |
2674 | | pos = line; |
2675 | | for (int i = 0; i < nAcross; i++) { |
2676 | | char * comma = strchr(pos, ','); |
2677 | | if (comma) { |
2678 | | *comma = '\0'; |
2679 | | } else if (i < nAcross - 1) { |
2680 | | nBadLine++; |
2681 | | break; |
2682 | | } |
2683 | | switch (orderColumn[i]) { |
2684 | | // name |
2685 | | case 0: |
2686 | | // For large problems this could be slow |
2687 | | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
2688 | | if (!strcmp(columnNames[iColumn], pos)) |
2689 | | break; |
2690 | | } |
2691 | | if (iColumn == numberColumns_) |
2692 | | iColumn = -1; |
2693 | | break; |
2694 | | // number |
2695 | | case 1: |
2696 | | iColumn = atoi(pos); |
2697 | | if (iColumn < 0 || iColumn >= numberColumns_) |
2698 | | iColumn = -1; |
2699 | | break; |
2700 | | // lower |
2701 | | case 2: |
2702 | | upper = atof(pos); |
2703 | | break; |
2704 | | // upper |
2705 | | case 3: |
2706 | | lower = atof(pos); |
2707 | | break; |
2708 | | // objective |
2709 | | case 4: |
2710 | | obj = atof(pos); |
2711 | | upper = lower; |
2712 | | break; |
2713 | | } |
2714 | | if (comma) { |
2715 | | *comma = ','; |
2716 | | pos = comma + 1; |
2717 | | } |
2718 | | } |
2719 | | if (iColumn >= 0) { |
2720 | | if (columnLower_[iColumn]>-1.0e20) |
2721 | | lowerColumnMove[iColumn] = lower; |
2722 | | else |
2723 | | lowerColumnMove[iColumn]=0.0; |
2724 | | if (columnUpper_[iColumn]<1.0e20) |
2725 | | upperColumnMove[iColumn] = upper; |
2726 | | else |
2727 | | upperColumnMove[iColumn] = lower; |
2728 | | objectiveMove[iColumn] = obj; |
2729 | | } else { |
2730 | | nBadName++; |
2731 | | if(saveLine[0]=='\0') |
2732 | | strcpy(saveLine,line); |
2733 | | } |
2734 | | } |
2735 | | sprintf(line,"%d Column fields and %d records", nAcross, nLine); |
2736 | | handler_->message(CLP_GENERAL,messages_) |
2737 | | << line << CoinMessageEol; |
2738 | | if (nBadName) { |
2739 | | sprintf(line," ** %d records did not match on name/sequence, first bad %s", nBadName,saveLine); |
2740 | | handler_->message(CLP_GENERAL,messages_) |
2741 | | << line << CoinMessageEol; |
2742 | | returnCode=-1; |
2743 | | good=false; |
2744 | | } |
2745 | | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
2746 | | free(columnNames[iColumn]); |
2747 | | } |
2748 | | delete [] columnNames; |
| 2612 | char **columnNames = new char *[numberColumns_]; |
| 2613 | int iColumn; |
| 2614 | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 2615 | columnNames[iColumn] = CoinStrdup(columnName(iColumn).c_str()); |
| 2616 | } |
| 2617 | lowerColumnMove = new double[numberColumns_]; |
| 2618 | memset(lowerColumnMove, 0, numberColumns_ * sizeof(double)); |
| 2619 | upperColumnMove = new double[numberColumns_]; |
| 2620 | memset(upperColumnMove, 0, numberColumns_ * sizeof(double)); |
| 2621 | objectiveMove = new double[numberColumns_]; |
| 2622 | memset(objectiveMove, 0, numberColumns_ * sizeof(double)); |
| 2623 | int nLine = 0; |
| 2624 | int nBadLine = 0; |
| 2625 | int nBadName = 0; |
| 2626 | while (fgets(line, 200, fp)) { |
| 2627 | if (!strncmp(line, "ENDATA", 6)) |
| 2628 | break; |
| 2629 | nLine++; |
| 2630 | iColumn = -1; |
| 2631 | double upper = 0.0; |
| 2632 | double lower = 0.0; |
| 2633 | double obj = 0.0; |
| 2634 | char *pos = line; |
| 2635 | char *put = line; |
| 2636 | while (*pos >= ' ' && *pos != '\n') { |
| 2637 | if (*pos != ' ' && *pos != '\t') { |
| 2638 | *put = *pos; |
| 2639 | put++; |
| 2640 | } |
| 2641 | pos++; |
| 2642 | } |
| 2643 | *put = '\0'; |
| 2644 | pos = line; |
| 2645 | for (int i = 0; i < nAcross; i++) { |
| 2646 | char *comma = strchr(pos, ','); |
| 2647 | if (comma) { |
| 2648 | *comma = '\0'; |
| 2649 | } else if (i < nAcross - 1) { |
| 2650 | nBadLine++; |
| 2651 | break; |
| 2652 | } |
| 2653 | switch (orderColumn[i]) { |
| 2654 | // name |
| 2655 | case 0: |
| 2656 | // For large problems this could be slow |
| 2657 | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 2658 | if (!strcmp(columnNames[iColumn], pos)) |
| 2659 | break; |
| 2660 | } |
| 2661 | if (iColumn == numberColumns_) |
| 2662 | iColumn = -1; |
| 2663 | break; |
| 2664 | // number |
| 2665 | case 1: |
| 2666 | iColumn = atoi(pos); |
| 2667 | if (iColumn < 0 || iColumn >= numberColumns_) |
| 2668 | iColumn = -1; |
| 2669 | break; |
| 2670 | // lower |
| 2671 | case 2: |
| 2672 | upper = atof(pos); |
| 2673 | break; |
| 2674 | // upper |
| 2675 | case 3: |
| 2676 | lower = atof(pos); |
| 2677 | break; |
| 2678 | // objective |
| 2679 | case 4: |
| 2680 | obj = atof(pos); |
| 2681 | upper = lower; |
| 2682 | break; |
| 2683 | } |
| 2684 | if (comma) { |
| 2685 | *comma = ','; |
| 2686 | pos = comma + 1; |
| 2687 | } |
| 2688 | } |
| 2689 | if (iColumn >= 0) { |
| 2690 | if (columnLower_[iColumn] > -1.0e20) |
| 2691 | lowerColumnMove[iColumn] = lower; |
| 2692 | else |
| 2693 | lowerColumnMove[iColumn] = 0.0; |
| 2694 | if (columnUpper_[iColumn] < 1.0e20) |
| 2695 | upperColumnMove[iColumn] = upper; |
| 2696 | else |
| 2697 | upperColumnMove[iColumn] = lower; |
| 2698 | objectiveMove[iColumn] = obj; |
| 2699 | } else { |
| 2700 | nBadName++; |
| 2701 | if (saveLine[0] == '\0') |
| 2702 | strcpy(saveLine, line); |
| 2703 | } |
| 2704 | } |
| 2705 | sprintf(line, "%d Column fields and %d records", nAcross, nLine); |
| 2706 | handler_->message(CLP_GENERAL, messages_) |
| 2707 | << line << CoinMessageEol; |
| 2708 | if (nBadName) { |
| 2709 | sprintf(line, " ** %d records did not match on name/sequence, first bad %s", nBadName, saveLine); |
| 2710 | handler_->message(CLP_GENERAL, messages_) |
| 2711 | << line << CoinMessageEol; |
| 2712 | returnCode = -1; |
| 2713 | good = false; |
| 2714 | } |
| 2715 | for (iColumn = 0; iColumn < numberColumns_; iColumn++) { |
| 2716 | free(columnNames[iColumn]); |
| 2717 | } |
| 2718 | delete[] columnNames; |
3069 | | memcpy(saveLower+unscaledChangesOffset,lowerChange,numberTotal*sizeof(double)); |
3070 | | memcpy(saveUpper+unscaledChangesOffset,upperChange,numberTotal*sizeof(double)); |
3071 | | int nLowerChange=0; |
3072 | | int nUpperChange=0; |
3073 | | for (int i=0;i<numberColumns_;i++) { |
3074 | | if (lowerChange[i]) { |
3075 | | lowerList[nLowerChange++]=i; |
3076 | | } |
3077 | | if (upperChange[i]) { |
3078 | | upperList[nUpperChange++]=i; |
3079 | | } |
3080 | | } |
3081 | | lowerList[-2]=nLowerChange; |
3082 | | upperList[-2]=nUpperChange; |
3083 | | for (int i=numberColumns_;i<numberTotal;i++) { |
3084 | | if (lowerChange[i]) { |
3085 | | lowerList[nLowerChange++]=i; |
3086 | | } |
3087 | | if (upperChange[i]) { |
3088 | | upperList[nUpperChange++]=i; |
3089 | | } |
3090 | | } |
3091 | | lowerList[-1]=nLowerChange; |
3092 | | upperList[-1]=nUpperChange; |
3093 | | memcpy(lowerCopy,columnLower_,numberColumns_*sizeof(double)); |
3094 | | memcpy(upperCopy,columnUpper_,numberColumns_*sizeof(double)); |
3095 | | memcpy(lowerCopy+numberColumns_, |
3096 | | rowLower_,numberRows_*sizeof(double)); |
3097 | | memcpy(upperCopy+numberColumns_, |
3098 | | rowUpper_,numberRows_*sizeof(double)); |
| 3036 | memcpy(saveLower + unscaledChangesOffset, lowerChange, numberTotal * sizeof(double)); |
| 3037 | memcpy(saveUpper + unscaledChangesOffset, upperChange, numberTotal * sizeof(double)); |
| 3038 | int nLowerChange = 0; |
| 3039 | int nUpperChange = 0; |
| 3040 | for (int i = 0; i < numberColumns_; i++) { |
| 3041 | if (lowerChange[i]) { |
| 3042 | lowerList[nLowerChange++] = i; |
| 3043 | } |
| 3044 | if (upperChange[i]) { |
| 3045 | upperList[nUpperChange++] = i; |
| 3046 | } |
| 3047 | } |
| 3048 | lowerList[-2] = nLowerChange; |
| 3049 | upperList[-2] = nUpperChange; |
| 3050 | for (int i = numberColumns_; i < numberTotal; i++) { |
| 3051 | if (lowerChange[i]) { |
| 3052 | lowerList[nLowerChange++] = i; |
| 3053 | } |
| 3054 | if (upperChange[i]) { |
| 3055 | upperList[nUpperChange++] = i; |
| 3056 | } |
| 3057 | } |
| 3058 | lowerList[-1] = nLowerChange; |
| 3059 | upperList[-1] = nUpperChange; |
| 3060 | memcpy(lowerCopy, columnLower_, numberColumns_ * sizeof(double)); |
| 3061 | memcpy(upperCopy, columnUpper_, numberColumns_ * sizeof(double)); |
| 3062 | memcpy(lowerCopy + numberColumns_, |
| 3063 | rowLower_, numberRows_ * sizeof(double)); |
| 3064 | memcpy(upperCopy + numberColumns_, |
| 3065 | rowUpper_, numberRows_ * sizeof(double)); |