From 6e5297e714ec9708c9a538e78966ef94333abe08 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Fri, 25 Oct 2024 13:32:41 +0100 Subject: [PATCH 1/8] Added presolve-slacks test case to TestPresolve.cpp --- check/TestPresolve.cpp | 38 +++++++++++++++++++++++++++++++++ src/qpsolver/dantzigpricing.hpp | 3 ++- src/qpsolver/devexpricing.hpp | 2 ++ 3 files changed, 42 insertions(+), 1 deletion(-) diff --git a/check/TestPresolve.cpp b/check/TestPresolve.cpp index 86c5dfc469..b0a4dbc95b 100644 --- a/check/TestPresolve.cpp +++ b/check/TestPresolve.cpp @@ -601,3 +601,41 @@ TEST_CASE("write-presolved-model", "[highs_test_presolve]") { REQUIRE(highs.getInfo().simplex_iteration_count == -1); std::remove(presolved_model_file.c_str()); } + +TEST_CASE("presolve-slacks", "[highs_test_presolve]") { + // This LP reduces to empty, because the equation is a doubleton + HighsLp lp; + lp.num_col_ = 2; + lp.num_row_ = 1; + lp.col_cost_ = {1, 0}; + lp.col_lower_ = {0, 0}; + lp.col_upper_ = {kHighsInf, kHighsInf}; + lp.row_lower_ = {1}; + lp.row_upper_ = {1}; + lp.a_matrix_.start_ = {0, 1, 2}; + lp.a_matrix_.index_ = {0, 0}; + lp.a_matrix_.value_ = {1, 1}; + Highs h; + h.setOptionValue("output_flag", dev_run); + REQUIRE(h.passModel(lp) == HighsStatus::kOk); + REQUIRE(h.presolve() == HighsStatus::kOk); + REQUIRE(h.getPresolvedLp().num_col_ == 0); + REQUIRE(h.getPresolvedLp().num_row_ == 0); + + lp.num_col_ = 4; + lp.num_row_ = 2; + lp.col_cost_ = {-10, -25, 0, 0}; + lp.col_lower_ = {0, 0, 0, 0}; + lp.col_upper_ = {kHighsInf, kHighsInf, kHighsInf, kHighsInf}; + lp.row_lower_ = {80, 120}; + lp.row_upper_ = {80, 120}; + lp.a_matrix_.start_ = {0, 2, 4, 5, 6}; + lp.a_matrix_.index_ = {0, 1, 0, 1, 0, 1}; + lp.a_matrix_.value_ = {1, 1, 2, 4, 1, 1}; + + REQUIRE(h.passModel(lp) == HighsStatus::kOk); + REQUIRE(h.run() == HighsStatus::kOk); + REQUIRE(h.presolve() == HighsStatus::kOk); + // REQUIRE(h.getPresolvedLp().num_col_ == 0); + // REQUIRE(h.getPresolvedLp().num_row_ == 0); +} diff --git a/src/qpsolver/dantzigpricing.hpp b/src/qpsolver/dantzigpricing.hpp index 89b109e5c1..5a62dce7fb 100644 --- a/src/qpsolver/dantzigpricing.hpp +++ b/src/qpsolver/dantzigpricing.hpp @@ -52,8 +52,9 @@ class DantzigPricing : public Pricing { public: DantzigPricing(Runtime& rt, Basis& bas, ReducedCosts& rc) + // clang-format off : runtime(rt), basis(bas), redcosts(rc) {}; - + // clang-format on HighsInt price(const QpVector& x, const QpVector& gradient) { HighsInt minidx = chooseconstrainttodrop(redcosts.getReducedCosts()); return minidx; diff --git a/src/qpsolver/devexpricing.hpp b/src/qpsolver/devexpricing.hpp index a88f3ecea5..014d8ce283 100644 --- a/src/qpsolver/devexpricing.hpp +++ b/src/qpsolver/devexpricing.hpp @@ -58,7 +58,9 @@ class DevexPricing : public Pricing { : runtime(rt), basis(bas), redcosts(rc), + // clang-format off weights(std::vector(rt.instance.num_var, 1.0)) {}; + // clang-format on // B lambda = g // lambda = inv(B)g From fe88b22201e1fef5a88ded9b11d29865db6163c6 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Fri, 8 Nov 2024 13:42:57 +0000 Subject: [PATCH 2/8] Added skeleton SlackColSubstitution as presolve/postsolve action --- check/TestPresolve.cpp | 6 ++-- src/lp_data/Highs.cpp | 36 +++++++++++++++++++- src/lp_data/HighsOptions.h | 7 ++++ src/presolve/HighsPostsolveStack.cpp | 48 ++++++++++++++++++++++++++ src/presolve/HighsPostsolveStack.h | 50 ++++++++++++++++++++++++++++ 5 files changed, 143 insertions(+), 4 deletions(-) diff --git a/check/TestPresolve.cpp b/check/TestPresolve.cpp index b0a4dbc95b..7b676dc011 100644 --- a/check/TestPresolve.cpp +++ b/check/TestPresolve.cpp @@ -616,7 +616,7 @@ TEST_CASE("presolve-slacks", "[highs_test_presolve]") { lp.a_matrix_.index_ = {0, 0}; lp.a_matrix_.value_ = {1, 1}; Highs h; - h.setOptionValue("output_flag", dev_run); + // h.setOptionValue("output_flag", dev_run); REQUIRE(h.passModel(lp) == HighsStatus::kOk); REQUIRE(h.presolve() == HighsStatus::kOk); REQUIRE(h.getPresolvedLp().num_col_ == 0); @@ -636,6 +636,6 @@ TEST_CASE("presolve-slacks", "[highs_test_presolve]") { REQUIRE(h.passModel(lp) == HighsStatus::kOk); REQUIRE(h.run() == HighsStatus::kOk); REQUIRE(h.presolve() == HighsStatus::kOk); - // REQUIRE(h.getPresolvedLp().num_col_ == 0); - // REQUIRE(h.getPresolvedLp().num_row_ == 0); + REQUIRE(h.getPresolvedLp().num_col_ == 2); + REQUIRE(h.getPresolvedLp().num_row_ == 2); } diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index c1f1f48b72..87d6309d6f 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -1251,6 +1251,40 @@ HighsStatus Highs::run() { // Run solver. bool have_optimal_solution = false; // ToDo Put solution of presolved problem in a separate method + + if (!this->options_.presolve_remove_slacks) { + HighsLp& reduced_lp = presolve_.getReducedProblem(); + HighsInt num_double_slack = 0; + HighsInt num_slack = 0; + HighsInt num_zero_cost_slack = 0; + HighsInt num_unit_coeff_slack = 0; + double min_slack_coeff = kHighsInf; + double max_slack_coeff = -kHighsInf; + std::vector found_slack; + found_slack.assign(reduced_lp.num_row_, false); + for (HighsInt iCol = 0; iCol < reduced_lp.num_col_; iCol++) { + HighsInt nnz = reduced_lp.a_matrix_.start_[iCol+1] - reduced_lp.a_matrix_.start_[iCol]; + if (nnz != 1) continue; + HighsInt iRow = reduced_lp.a_matrix_.index_[reduced_lp.a_matrix_.start_[iCol]]; + if (found_slack[iRow]) { + num_double_slack++; + continue; + } + if (reduced_lp.row_lower_[iRow] != reduced_lp.row_upper_[iRow]) continue; + double coeff = std::fabs(reduced_lp.a_matrix_.value_[reduced_lp.a_matrix_.start_[iCol]]); + if (coeff == 1.0) num_unit_coeff_slack++; + min_slack_coeff = std::min(coeff,min_slack_coeff); + max_slack_coeff = std::max(coeff,max_slack_coeff); + found_slack[iRow] = true; + num_slack++; + if (reduced_lp.col_cost_[iCol] == 0) num_zero_cost_slack++; + } + printf("grepSlack,model,col,slack,unit coeff,zero_cost,double,min coeff, max_coeff\n"); + printf("grepSlack,%s,%d, %d, %d, %d, %d, %g, %g\n", this->model_.lp_.model_name_.c_str(), + int(reduced_lp.num_col_), int(num_slack), int(num_unit_coeff_slack), int(num_zero_cost_slack), int(num_double_slack), + min_slack_coeff, max_slack_coeff); + } + switch (model_presolve_status_) { case HighsPresolveStatus::kNotPresolved: { ekk_instance_.lp_name_ = "Original LP"; @@ -1275,7 +1309,7 @@ HighsStatus Highs::run() { break; } case HighsPresolveStatus::kReduced: { - HighsLp& reduced_lp = presolve_.getReducedProblem(); + HighsLp& reduced_lp = presolve_.getReducedProblem(); reduced_lp.setMatrixDimensions(); if (kAllowDeveloperAssert) { // Validate the reduced LP diff --git a/src/lp_data/HighsOptions.h b/src/lp_data/HighsOptions.h index 6e2be6cc78..8103f4db7b 100644 --- a/src/lp_data/HighsOptions.h +++ b/src/lp_data/HighsOptions.h @@ -372,6 +372,7 @@ struct HighsOptionsStruct { HighsInt presolve_substitution_maxfillin; HighsInt presolve_rule_off; bool presolve_rule_logging; + bool presolve_remove_slacks; bool simplex_initial_condition_check; bool no_unnecessary_rebuild_refactor; double simplex_initial_condition_tolerance; @@ -507,6 +508,7 @@ struct HighsOptionsStruct { presolve_substitution_maxfillin(0), presolve_rule_off(0), presolve_rule_logging(false), + presolve_remove_slacks(false), simplex_initial_condition_check(false), no_unnecessary_rebuild_refactor(false), simplex_initial_condition_tolerance(0.0), @@ -1324,6 +1326,11 @@ class HighsOptions : public HighsOptionsStruct { advanced, &presolve_rule_logging, false); records.push_back(record_bool); + record_bool = new OptionRecordBool( + "presolve_remove_slacks", "Remove slacks after presolve", + advanced, &presolve_remove_slacks, true);//false); + records.push_back(record_bool); + record_int = new OptionRecordInt( "presolve_substitution_maxfillin", "Maximal fillin allowed for substitutions in presolve", advanced, diff --git a/src/presolve/HighsPostsolveStack.cpp b/src/presolve/HighsPostsolveStack.cpp index bec3e4945f..79d582a7f4 100644 --- a/src/presolve/HighsPostsolveStack.cpp +++ b/src/presolve/HighsPostsolveStack.cpp @@ -1351,4 +1351,52 @@ void HighsPostsolveStack::DuplicateColumn::transformToPresolvedSpace( primalSol[col] = primalSol[col] + colScale * primalSol[duplicateCol]; } +void HighsPostsolveStack::SlackColSubstitution::undo( + const HighsOptions& options, const std::vector& rowValues, + const std::vector& colValues, HighsSolution& solution, + HighsBasis& basis) { + // a (removed) cut may have been used in this reduction. + bool isModelRow = static_cast(row) < solution.row_value.size(); + + // compute primal values + double colCoef = 0; + HighsCDouble rowValue = 0; + for (const auto& rowVal : rowValues) { + if (rowVal.index == col) + colCoef = rowVal.value; + else + rowValue += rowVal.value * solution.col_value[rowVal.index]; + } + + assert(colCoef != 0); + // Row values aren't fully postsolved, so why do this? + if (isModelRow) + solution.row_value[row] = + double(rowValue + colCoef * solution.col_value[col]); + solution.col_value[col] = double((rhs - rowValue) / colCoef); + + // if no dual values requested, return here + if (!solution.dual_valid) return; + + // compute the row dual value such that reduced cost of basic column is 0 + if (isModelRow) { + solution.row_dual[row] = 0; + HighsCDouble dualval = colCost; + for (const auto& colVal : colValues) { + if (static_cast(colVal.index) < solution.row_dual.size()) + dualval -= colVal.value * solution.row_dual[colVal.index]; + } + solution.row_dual[row] = double(dualval / colCoef); + } + + solution.col_dual[col] = 0; + + // set basis status if necessary + if (!basis.valid) return; + + basis.col_status[col] = HighsBasisStatus::kBasic; + if (isModelRow) + basis.row_status[row] = computeRowStatus(solution.row_dual[row], rowType); +} + } // namespace presolve diff --git a/src/presolve/HighsPostsolveStack.h b/src/presolve/HighsPostsolveStack.h index 2200977b5b..79a658f977 100644 --- a/src/presolve/HighsPostsolveStack.h +++ b/src/presolve/HighsPostsolveStack.h @@ -219,6 +219,19 @@ class HighsPostsolveStack { void transformToPresolvedSpace(std::vector& primalSol) const; }; + struct SlackColSubstitution { + double rhs; + double colCost; + HighsInt row; + HighsInt col; + RowType rowType; + + void undo(const HighsOptions& options, + const std::vector& rowValues, + const std::vector& colValues, HighsSolution& solution, + HighsBasis& basis); + }; + /// tags for reduction enum class ReductionType : uint8_t { kLinearTransform, @@ -234,6 +247,7 @@ class HighsPostsolveStack { kForcingColumnRemovedRow, kDuplicateRow, kDuplicateColumn, + kSlackColSubstitution, }; HighsDataStack reductionValues; @@ -323,6 +337,26 @@ class HighsPostsolveStack { reductionAdded(ReductionType::kFreeColSubstitution); } + template + void slackColSubstitution(HighsInt row, HighsInt col, double rhs, + double colCost, RowType rowType, + const HighsMatrixSlice& rowVec, + const HighsMatrixSlice& colVec) { + rowValues.clear(); + for (const HighsSliceNonzero& rowVal : rowVec) + rowValues.emplace_back(origColIndex[rowVal.index()], rowVal.value()); + + colValues.clear(); + for (const HighsSliceNonzero& colVal : colVec) + colValues.emplace_back(origRowIndex[colVal.index()], colVal.value()); + + reductionValues.push(SlackColSubstitution{rhs, colCost, origRowIndex[row], + origColIndex[col], rowType}); + reductionValues.push(rowValues); + reductionValues.push(colValues); + reductionAdded(ReductionType::kSlackColSubstitution); + } + template void doubletonEquation(HighsInt row, HighsInt colSubst, HighsInt col, double coefSubst, double coef, double rhs, @@ -710,6 +744,14 @@ class HighsPostsolveStack { reduction.undo(options, solution, basis); break; } + case ReductionType::kSlackColSubstitution: { + SlackColSubstitution reduction; + reductionValues.pop(colValues); + reductionValues.pop(rowValues); + reductionValues.pop(reduction); + reduction.undo(options, rowValues, colValues, solution, basis); + break; + } default: printf("Reduction case %d not handled\n", int(reductions[i - 1].first)); @@ -887,6 +929,14 @@ class HighsPostsolveStack { reductionValues.pop(reduction); reduction.undo(options, solution, basis); } + case ReductionType::kSlackColSubstitution: { + SlackColSubstitution reduction; + reductionValues.pop(colValues); + reductionValues.pop(rowValues); + reductionValues.pop(reduction); + reduction.undo(options, rowValues, colValues, solution, basis); + break; + } } } #ifdef DEBUG_EXTRA From 14c13c979d66dac55a5745361c1594c0b63e180a Mon Sep 17 00:00:00 2001 From: JAJHall Date: Fri, 8 Nov 2024 13:44:13 +0000 Subject: [PATCH 3/8] Switched off presolve-slacks check for elimination of slacks --- check/TestPresolve.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/check/TestPresolve.cpp b/check/TestPresolve.cpp index 7b676dc011..40d1ea900c 100644 --- a/check/TestPresolve.cpp +++ b/check/TestPresolve.cpp @@ -636,6 +636,6 @@ TEST_CASE("presolve-slacks", "[highs_test_presolve]") { REQUIRE(h.passModel(lp) == HighsStatus::kOk); REQUIRE(h.run() == HighsStatus::kOk); REQUIRE(h.presolve() == HighsStatus::kOk); - REQUIRE(h.getPresolvedLp().num_col_ == 2); - REQUIRE(h.getPresolvedLp().num_row_ == 2); + // REQUIRE(h.getPresolvedLp().num_col_ == 2); + // REQUIRE(h.getPresolvedLp().num_row_ == 2); } From d5ffacc534dc52201b1aee0ed87f5e819814cd48 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Fri, 8 Nov 2024 15:46:11 +0000 Subject: [PATCH 4/8] Now getting through to postsolve after eliminating a slack in presolve for adlittle --- src/lp_data/Highs.cpp | 7 ++-- src/presolve/HPresolve.cpp | 61 ++++++++++++++++++++++++++++ src/presolve/HPresolve.h | 2 + src/presolve/HighsPostsolveStack.cpp | 7 +++- src/presolve/HighsPostsolveStack.h | 19 ++++----- 5 files changed, 81 insertions(+), 15 deletions(-) diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index 87d6309d6f..05d89e3ca4 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -1252,7 +1252,7 @@ HighsStatus Highs::run() { bool have_optimal_solution = false; // ToDo Put solution of presolved problem in a separate method - if (!this->options_.presolve_remove_slacks) { + // if (!this->options_.presolve_remove_slacks) { HighsLp& reduced_lp = presolve_.getReducedProblem(); HighsInt num_double_slack = 0; HighsInt num_slack = 0; @@ -1271,19 +1271,20 @@ HighsStatus Highs::run() { continue; } if (reduced_lp.row_lower_[iRow] != reduced_lp.row_upper_[iRow]) continue; + num_slack++; + printf("Column %d is slack\n", int(iCol)); double coeff = std::fabs(reduced_lp.a_matrix_.value_[reduced_lp.a_matrix_.start_[iCol]]); if (coeff == 1.0) num_unit_coeff_slack++; min_slack_coeff = std::min(coeff,min_slack_coeff); max_slack_coeff = std::max(coeff,max_slack_coeff); found_slack[iRow] = true; - num_slack++; if (reduced_lp.col_cost_[iCol] == 0) num_zero_cost_slack++; } printf("grepSlack,model,col,slack,unit coeff,zero_cost,double,min coeff, max_coeff\n"); printf("grepSlack,%s,%d, %d, %d, %d, %d, %g, %g\n", this->model_.lp_.model_name_.c_str(), int(reduced_lp.num_col_), int(num_slack), int(num_unit_coeff_slack), int(num_zero_cost_slack), int(num_double_slack), min_slack_coeff, max_slack_coeff); - } + // } switch (model_presolve_status_) { case HighsPresolveStatus::kNotPresolved: { diff --git a/src/presolve/HPresolve.cpp b/src/presolve/HPresolve.cpp index 5551ef0f20..18f6df36c8 100644 --- a/src/presolve/HPresolve.cpp +++ b/src/presolve/HPresolve.cpp @@ -4360,6 +4360,11 @@ HPresolve::Result HPresolve::presolve(HighsPostsolveStack& postsolve_stack) { break; } + // Now consider removing slacks + if (options->presolve_remove_slacks) { + HPRESOLVE_CHECKED_CALL(removeSlacks(postsolve_stack)); + } + report(); } else { highsLogUser(options->log_options, HighsLogType::kInfo, @@ -4375,6 +4380,62 @@ HPresolve::Result HPresolve::presolve(HighsPostsolveStack& postsolve_stack) { return Result::kOk; } +HPresolve::Result HPresolve::removeSlacks(HighsPostsolveStack& postsolve_stack) { + // singletonColumns data structure appears not to be retained + // throughout presolve + // + bool unit_coeff_only = true; + for (HighsInt iCol = 0; iCol != model->num_col_; ++iCol) { + if (colDeleted[iCol]) continue; + if (colsize[iCol] != 1) continue; + if (model->integrality_[iCol] == HighsVarType::kInteger) continue; + HighsInt coliter = colhead[iCol]; + HighsInt iRow = Arow[coliter]; + assert(Acol[coliter] == iCol); + assert(!rowDeleted[iRow]); + if (model->row_lower_[iRow] != model->row_upper_[iRow]) continue; + double lower = model->col_lower_[iCol]; + double upper = model->col_upper_[iCol]; + double cost = model->col_cost_[iCol]; + double rhs = model->row_lower_[iRow]; + double coeff = Avalue[coliter]; + printf("Col %d is continuous and is singleton in equality row %d with cost %g, bounds [%g, %g], coeff %g and RHS = %g\n", int(iCol), int(iRow), cost, lower, upper, coeff, rhs); + if (unit_coeff_only && std::fabs(coeff) != 1.0) continue; + assert(coeff); + // Slack is s = (rhs - a^Tx)/coeff + // + if (coeff > 0) { + // Constraint bounds become [rhs - coeff * lower, rhs - coeff * + // upper] + model->col_lower_[iCol] = rhs - coeff * lower; + model->col_upper_[iCol] = rhs - coeff * upper; + } else { + // Constraint bounds become [rhs - coeff * upper, rhs - coeff * + // lower] + model->col_lower_[iCol] = rhs - coeff * upper; + model->col_upper_[iCol] = rhs - coeff * lower; + } + if (cost) { + // Cost is (cost * rhs / coeff) + (col_cost - (cost/coeff) row_values)^Tx + double multiplier = cost / coeff; + for (const HighsSliceNonzero& nonzero : getRowVector(iRow)) { + HighsInt local_iCol = nonzero.index(); + double local_value = nonzero.value(); + model->col_cost_[local_iCol] -= multiplier * local_value; + } + model->offset_ += multiplier * rhs; + } + // + postsolve_stack.slackColSubstitution(iRow, iCol, rhs, cost, lower, upper, coeff, + getColumnVector(iCol)); + markColDeleted(iCol); + + unlink(coliter); + + } + return Result::kOk; +} + HPresolve::Result HPresolve::checkLimits(HighsPostsolveStack& postsolve_stack) { size_t numreductions = postsolve_stack.numReductions(); diff --git a/src/presolve/HPresolve.h b/src/presolve/HPresolve.h index 2a5b3427c8..28ac2bbe47 100644 --- a/src/presolve/HPresolve.h +++ b/src/presolve/HPresolve.h @@ -269,6 +269,8 @@ class HPresolve { Result presolve(HighsPostsolveStack& postsolve_stack); + Result removeSlacks(HighsPostsolveStack& postsolve_stack); + Result checkLimits(HighsPostsolveStack& postsolve_stack); void storeCurrentProblemSize(); diff --git a/src/presolve/HighsPostsolveStack.cpp b/src/presolve/HighsPostsolveStack.cpp index 79d582a7f4..a2bc80f2e2 100644 --- a/src/presolve/HighsPostsolveStack.cpp +++ b/src/presolve/HighsPostsolveStack.cpp @@ -1355,6 +1355,11 @@ void HighsPostsolveStack::SlackColSubstitution::undo( const HighsOptions& options, const std::vector& rowValues, const std::vector& colValues, HighsSolution& solution, HighsBasis& basis) { + + assert(111==222); + + + // a (removed) cut may have been used in this reduction. bool isModelRow = static_cast(row) < solution.row_value.size(); @@ -1396,7 +1401,7 @@ void HighsPostsolveStack::SlackColSubstitution::undo( basis.col_status[col] = HighsBasisStatus::kBasic; if (isModelRow) - basis.row_status[row] = computeRowStatus(solution.row_dual[row], rowType); + basis.row_status[row] = computeRowStatus(solution.row_dual[row], RowType::kEq); } } // namespace presolve diff --git a/src/presolve/HighsPostsolveStack.h b/src/presolve/HighsPostsolveStack.h index 79a658f977..c4160dd226 100644 --- a/src/presolve/HighsPostsolveStack.h +++ b/src/presolve/HighsPostsolveStack.h @@ -222,9 +222,11 @@ class HighsPostsolveStack { struct SlackColSubstitution { double rhs; double colCost; + double colLower; + double colUpper; + double colCoeff; HighsInt row; HighsInt col; - RowType rowType; void undo(const HighsOptions& options, const std::vector& rowValues, @@ -337,21 +339,16 @@ class HighsPostsolveStack { reductionAdded(ReductionType::kFreeColSubstitution); } - template + template void slackColSubstitution(HighsInt row, HighsInt col, double rhs, - double colCost, RowType rowType, - const HighsMatrixSlice& rowVec, - const HighsMatrixSlice& colVec) { - rowValues.clear(); - for (const HighsSliceNonzero& rowVal : rowVec) - rowValues.emplace_back(origColIndex[rowVal.index()], rowVal.value()); - + double colCost, double colLower, double colUpper, double colCoeff, + const HighsMatrixSlice& colVec) { colValues.clear(); for (const HighsSliceNonzero& colVal : colVec) colValues.emplace_back(origRowIndex[colVal.index()], colVal.value()); - reductionValues.push(SlackColSubstitution{rhs, colCost, origRowIndex[row], - origColIndex[col], rowType}); + reductionValues.push(SlackColSubstitution{rhs, colCost, colLower, colUpper, colCoeff, origRowIndex[row], + origColIndex[col]}); reductionValues.push(rowValues); reductionValues.push(colValues); reductionAdded(ReductionType::kSlackColSubstitution); From b14bbe1a49fb08bc10dc4a7dc86b93792650e15d Mon Sep 17 00:00:00 2001 From: JAJHall Date: Fri, 8 Nov 2024 17:43:16 +0000 Subject: [PATCH 5/8] Solves adlittle with slack removed --- src/presolve/HPresolve.cpp | 23 ++++++++++------------- src/presolve/HighsPostsolveStack.cpp | 15 +++++++++------ src/presolve/HighsPostsolveStack.h | 16 +++++++++++----- 3 files changed, 30 insertions(+), 24 deletions(-) diff --git a/src/presolve/HPresolve.cpp b/src/presolve/HPresolve.cpp index 18f6df36c8..e1e721347c 100644 --- a/src/presolve/HPresolve.cpp +++ b/src/presolve/HPresolve.cpp @@ -4404,17 +4404,13 @@ HPresolve::Result HPresolve::removeSlacks(HighsPostsolveStack& postsolve_stack) assert(coeff); // Slack is s = (rhs - a^Tx)/coeff // - if (coeff > 0) { - // Constraint bounds become [rhs - coeff * lower, rhs - coeff * - // upper] - model->col_lower_[iCol] = rhs - coeff * lower; - model->col_upper_[iCol] = rhs - coeff * upper; - } else { - // Constraint bounds become [rhs - coeff * upper, rhs - coeff * - // lower] - model->col_lower_[iCol] = rhs - coeff * upper; - model->col_upper_[iCol] = rhs - coeff * lower; - } + // Constraint bounds become: + // + // For coeff > 0 [rhs - coeff * upper, rhs - coeff * lower] + // + // For coeff < 0 [rhs - coeff * lower, rhs - coeff * upper] + model->row_lower_[iRow] = coeff > 0 ? rhs - coeff * upper : rhs - coeff * lower; + model->row_upper_[iRow] = coeff > 0 ? rhs - coeff * lower : rhs - coeff * upper; if (cost) { // Cost is (cost * rhs / coeff) + (col_cost - (cost/coeff) row_values)^Tx double multiplier = cost / coeff; @@ -4426,8 +4422,9 @@ HPresolve::Result HPresolve::removeSlacks(HighsPostsolveStack& postsolve_stack) model->offset_ += multiplier * rhs; } // - postsolve_stack.slackColSubstitution(iRow, iCol, rhs, cost, lower, upper, coeff, - getColumnVector(iCol)); + postsolve_stack.slackColSubstitution(iRow, iCol, rhs, cost, lower, upper, //coeff, + getRowVector(iRow), + getColumnVector(iCol)); markColDeleted(iCol); unlink(coliter); diff --git a/src/presolve/HighsPostsolveStack.cpp b/src/presolve/HighsPostsolveStack.cpp index a2bc80f2e2..9a976187f9 100644 --- a/src/presolve/HighsPostsolveStack.cpp +++ b/src/presolve/HighsPostsolveStack.cpp @@ -1356,11 +1356,14 @@ void HighsPostsolveStack::SlackColSubstitution::undo( const std::vector& colValues, HighsSolution& solution, HighsBasis& basis) { - assert(111==222); - - - + // Taken from HighsPostsolveStack::FreeColSubstitution::undo( + // // a (removed) cut may have been used in this reduction. + // + // May have to determine row dual and basis status unless doing + // primal-only transformation in MIP solver, in which case row may + // no longer exist if it corresponds to a removed cut, so have to + // avoid exceeding array bounds of solution.row_value bool isModelRow = static_cast(row) < solution.row_value.size(); // compute primal values @@ -1375,9 +1378,9 @@ void HighsPostsolveStack::SlackColSubstitution::undo( assert(colCoef != 0); // Row values aren't fully postsolved, so why do this? - if (isModelRow) - solution.row_value[row] = + if (isModelRow) solution.row_value[row] = double(rowValue + colCoef * solution.col_value[col]); + printf("HighsPostsolveStack::SlackColSubstitution::undo rowValue = %g\n", double(rowValue)); solution.col_value[col] = double((rhs - rowValue) / colCoef); // if no dual values requested, return here diff --git a/src/presolve/HighsPostsolveStack.h b/src/presolve/HighsPostsolveStack.h index c4160dd226..954ffd2732 100644 --- a/src/presolve/HighsPostsolveStack.h +++ b/src/presolve/HighsPostsolveStack.h @@ -224,7 +224,7 @@ class HighsPostsolveStack { double colCost; double colLower; double colUpper; - double colCoeff; + // double colCoeff; HighsInt row; HighsInt col; @@ -339,16 +339,22 @@ class HighsPostsolveStack { reductionAdded(ReductionType::kFreeColSubstitution); } - template + template void slackColSubstitution(HighsInt row, HighsInt col, double rhs, - double colCost, double colLower, double colUpper, double colCoeff, + double colCost, double colLower, double colUpper, //double colCoeff, + const HighsMatrixSlice& rowVec, const HighsMatrixSlice& colVec) { + rowValues.clear(); + for (const HighsSliceNonzero& rowVal : rowVec) + rowValues.emplace_back(origColIndex[rowVal.index()], rowVal.value()); + colValues.clear(); for (const HighsSliceNonzero& colVal : colVec) colValues.emplace_back(origRowIndex[colVal.index()], colVal.value()); - reductionValues.push(SlackColSubstitution{rhs, colCost, colLower, colUpper, colCoeff, origRowIndex[row], - origColIndex[col]}); + reductionValues.push(SlackColSubstitution{rhs, colCost, colLower, colUpper, //colCoeff, + origRowIndex[row], + origColIndex[col]}); reductionValues.push(rowValues); reductionValues.push(colValues); reductionAdded(ReductionType::kSlackColSubstitution); From b32cfd0892d48b55ff58afbd799e28b052aff118 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Fri, 8 Nov 2024 18:21:37 +0000 Subject: [PATCH 6/8] Tidied and formatted --- check/TestPresolve.cpp | 8 ++-- src/lp_data/Highs.cpp | 72 +++++++++++++++------------- src/lp_data/HighsOptions.h | 8 ++-- src/presolve/HPresolve.cpp | 29 ++++++----- src/presolve/HPresolve.h | 2 +- src/presolve/HighsPostsolveStack.cpp | 34 ++++++++----- src/presolve/HighsPostsolveStack.h | 29 ++++------- src/qpsolver/dantzigpricing.hpp | 2 +- src/qpsolver/devexpricing.hpp | 2 +- 9 files changed, 98 insertions(+), 88 deletions(-) diff --git a/check/TestPresolve.cpp b/check/TestPresolve.cpp index 40d1ea900c..7d2e5d9604 100644 --- a/check/TestPresolve.cpp +++ b/check/TestPresolve.cpp @@ -616,7 +616,7 @@ TEST_CASE("presolve-slacks", "[highs_test_presolve]") { lp.a_matrix_.index_ = {0, 0}; lp.a_matrix_.value_ = {1, 1}; Highs h; - // h.setOptionValue("output_flag", dev_run); + h.setOptionValue("output_flag", dev_run); REQUIRE(h.passModel(lp) == HighsStatus::kOk); REQUIRE(h.presolve() == HighsStatus::kOk); REQUIRE(h.getPresolvedLp().num_col_ == 0); @@ -632,10 +632,10 @@ TEST_CASE("presolve-slacks", "[highs_test_presolve]") { lp.a_matrix_.start_ = {0, 2, 4, 5, 6}; lp.a_matrix_.index_ = {0, 1, 0, 1, 0, 1}; lp.a_matrix_.value_ = {1, 1, 2, 4, 1, 1}; - + REQUIRE(h.setOptionValue("presolve_remove_slacks", true) == HighsStatus::kOk); REQUIRE(h.passModel(lp) == HighsStatus::kOk); REQUIRE(h.run() == HighsStatus::kOk); REQUIRE(h.presolve() == HighsStatus::kOk); - // REQUIRE(h.getPresolvedLp().num_col_ == 2); - // REQUIRE(h.getPresolvedLp().num_row_ == 2); + REQUIRE(h.getPresolvedLp().num_col_ == 2); + REQUIRE(h.getPresolvedLp().num_row_ == 2); } diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index 05d89e3ca4..894b703aac 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -1251,40 +1251,46 @@ HighsStatus Highs::run() { // Run solver. bool have_optimal_solution = false; // ToDo Put solution of presolved problem in a separate method - + // if (!this->options_.presolve_remove_slacks) { - HighsLp& reduced_lp = presolve_.getReducedProblem(); - HighsInt num_double_slack = 0; - HighsInt num_slack = 0; - HighsInt num_zero_cost_slack = 0; - HighsInt num_unit_coeff_slack = 0; - double min_slack_coeff = kHighsInf; - double max_slack_coeff = -kHighsInf; - std::vector found_slack; - found_slack.assign(reduced_lp.num_row_, false); - for (HighsInt iCol = 0; iCol < reduced_lp.num_col_; iCol++) { - HighsInt nnz = reduced_lp.a_matrix_.start_[iCol+1] - reduced_lp.a_matrix_.start_[iCol]; - if (nnz != 1) continue; - HighsInt iRow = reduced_lp.a_matrix_.index_[reduced_lp.a_matrix_.start_[iCol]]; - if (found_slack[iRow]) { - num_double_slack++; - continue; - } - if (reduced_lp.row_lower_[iRow] != reduced_lp.row_upper_[iRow]) continue; - num_slack++; - printf("Column %d is slack\n", int(iCol)); - double coeff = std::fabs(reduced_lp.a_matrix_.value_[reduced_lp.a_matrix_.start_[iCol]]); - if (coeff == 1.0) num_unit_coeff_slack++; - min_slack_coeff = std::min(coeff,min_slack_coeff); - max_slack_coeff = std::max(coeff,max_slack_coeff); - found_slack[iRow] = true; - if (reduced_lp.col_cost_[iCol] == 0) num_zero_cost_slack++; + HighsLp& reduced_lp = presolve_.getReducedProblem(); + HighsInt num_double_slack = 0; + HighsInt num_slack = 0; + HighsInt num_zero_cost_slack = 0; + HighsInt num_unit_coeff_slack = 0; + double min_slack_coeff = kHighsInf; + double max_slack_coeff = -kHighsInf; + std::vector found_slack; + found_slack.assign(reduced_lp.num_row_, false); + for (HighsInt iCol = 0; iCol < reduced_lp.num_col_; iCol++) { + HighsInt nnz = reduced_lp.a_matrix_.start_[iCol + 1] - + reduced_lp.a_matrix_.start_[iCol]; + if (nnz != 1) continue; + HighsInt iRow = + reduced_lp.a_matrix_.index_[reduced_lp.a_matrix_.start_[iCol]]; + if (found_slack[iRow]) { + num_double_slack++; + continue; } - printf("grepSlack,model,col,slack,unit coeff,zero_cost,double,min coeff, max_coeff\n"); - printf("grepSlack,%s,%d, %d, %d, %d, %d, %g, %g\n", this->model_.lp_.model_name_.c_str(), - int(reduced_lp.num_col_), int(num_slack), int(num_unit_coeff_slack), int(num_zero_cost_slack), int(num_double_slack), - min_slack_coeff, max_slack_coeff); - // } + if (reduced_lp.row_lower_[iRow] != reduced_lp.row_upper_[iRow]) continue; + num_slack++; + printf("Column %d is slack\n", int(iCol)); + double coeff = std::fabs( + reduced_lp.a_matrix_.value_[reduced_lp.a_matrix_.start_[iCol]]); + if (coeff == 1.0) num_unit_coeff_slack++; + min_slack_coeff = std::min(coeff, min_slack_coeff); + max_slack_coeff = std::max(coeff, max_slack_coeff); + found_slack[iRow] = true; + if (reduced_lp.col_cost_[iCol] == 0) num_zero_cost_slack++; + } + printf( + "grepSlack,model,col,slack,unit coeff,zero_cost,double,min coeff, " + "max_coeff\n"); + printf("grepSlack,%s,%d, %d, %d, %d, %d, %g, %g\n", + this->model_.lp_.model_name_.c_str(), int(reduced_lp.num_col_), + int(num_slack), int(num_unit_coeff_slack), int(num_zero_cost_slack), + int(num_double_slack), min_slack_coeff, max_slack_coeff); + // } switch (model_presolve_status_) { case HighsPresolveStatus::kNotPresolved: { @@ -1310,7 +1316,7 @@ HighsStatus Highs::run() { break; } case HighsPresolveStatus::kReduced: { - HighsLp& reduced_lp = presolve_.getReducedProblem(); + HighsLp& reduced_lp = presolve_.getReducedProblem(); reduced_lp.setMatrixDimensions(); if (kAllowDeveloperAssert) { // Validate the reduced LP diff --git a/src/lp_data/HighsOptions.h b/src/lp_data/HighsOptions.h index 8103f4db7b..64d565014e 100644 --- a/src/lp_data/HighsOptions.h +++ b/src/lp_data/HighsOptions.h @@ -508,7 +508,7 @@ struct HighsOptionsStruct { presolve_substitution_maxfillin(0), presolve_rule_off(0), presolve_rule_logging(false), - presolve_remove_slacks(false), + presolve_remove_slacks(false), simplex_initial_condition_check(false), no_unnecessary_rebuild_refactor(false), simplex_initial_condition_tolerance(0.0), @@ -1326,9 +1326,9 @@ class HighsOptions : public HighsOptionsStruct { advanced, &presolve_rule_logging, false); records.push_back(record_bool); - record_bool = new OptionRecordBool( - "presolve_remove_slacks", "Remove slacks after presolve", - advanced, &presolve_remove_slacks, true);//false); + record_bool = new OptionRecordBool("presolve_remove_slacks", + "Remove slacks after presolve", advanced, + &presolve_remove_slacks, false); records.push_back(record_bool); record_int = new OptionRecordInt( diff --git a/src/presolve/HPresolve.cpp b/src/presolve/HPresolve.cpp index e1e721347c..de5a26204e 100644 --- a/src/presolve/HPresolve.cpp +++ b/src/presolve/HPresolve.cpp @@ -4380,7 +4380,8 @@ HPresolve::Result HPresolve::presolve(HighsPostsolveStack& postsolve_stack) { return Result::kOk; } -HPresolve::Result HPresolve::removeSlacks(HighsPostsolveStack& postsolve_stack) { +HPresolve::Result HPresolve::removeSlacks( + HighsPostsolveStack& postsolve_stack) { // singletonColumns data structure appears not to be retained // throughout presolve // @@ -4399,7 +4400,10 @@ HPresolve::Result HPresolve::removeSlacks(HighsPostsolveStack& postsolve_stack) double cost = model->col_cost_[iCol]; double rhs = model->row_lower_[iRow]; double coeff = Avalue[coliter]; - printf("Col %d is continuous and is singleton in equality row %d with cost %g, bounds [%g, %g], coeff %g and RHS = %g\n", int(iCol), int(iRow), cost, lower, upper, coeff, rhs); + printf( + "Col %d is continuous and is singleton in equality row %d with cost " + "%g, bounds [%g, %g], coeff %g and RHS = %g\n", + int(iCol), int(iRow), cost, lower, upper, coeff, rhs); if (unit_coeff_only && std::fabs(coeff) != 1.0) continue; assert(coeff); // Slack is s = (rhs - a^Tx)/coeff @@ -4409,26 +4413,27 @@ HPresolve::Result HPresolve::removeSlacks(HighsPostsolveStack& postsolve_stack) // For coeff > 0 [rhs - coeff * upper, rhs - coeff * lower] // // For coeff < 0 [rhs - coeff * lower, rhs - coeff * upper] - model->row_lower_[iRow] = coeff > 0 ? rhs - coeff * upper : rhs - coeff * lower; - model->row_upper_[iRow] = coeff > 0 ? rhs - coeff * lower : rhs - coeff * upper; + model->row_lower_[iRow] = + coeff > 0 ? rhs - coeff * upper : rhs - coeff * lower; + model->row_upper_[iRow] = + coeff > 0 ? rhs - coeff * lower : rhs - coeff * upper; if (cost) { // Cost is (cost * rhs / coeff) + (col_cost - (cost/coeff) row_values)^Tx double multiplier = cost / coeff; for (const HighsSliceNonzero& nonzero : getRowVector(iRow)) { - HighsInt local_iCol = nonzero.index(); - double local_value = nonzero.value(); - model->col_cost_[local_iCol] -= multiplier * local_value; + HighsInt local_iCol = nonzero.index(); + double local_value = nonzero.value(); + model->col_cost_[local_iCol] -= multiplier * local_value; } model->offset_ += multiplier * rhs; } - // - postsolve_stack.slackColSubstitution(iRow, iCol, rhs, cost, lower, upper, //coeff, - getRowVector(iRow), - getColumnVector(iCol)); + // + postsolve_stack.slackColSubstitution(iRow, iCol, rhs, cost, + getRowVector(iRow)); + markColDeleted(iCol); unlink(coliter); - } return Result::kOk; } diff --git a/src/presolve/HPresolve.h b/src/presolve/HPresolve.h index 28ac2bbe47..3699d23399 100644 --- a/src/presolve/HPresolve.h +++ b/src/presolve/HPresolve.h @@ -270,7 +270,7 @@ class HPresolve { Result presolve(HighsPostsolveStack& postsolve_stack); Result removeSlacks(HighsPostsolveStack& postsolve_stack); - + Result checkLimits(HighsPostsolveStack& postsolve_stack); void storeCurrentProblemSize(); diff --git a/src/presolve/HighsPostsolveStack.cpp b/src/presolve/HighsPostsolveStack.cpp index 9a976187f9..43485006ea 100644 --- a/src/presolve/HighsPostsolveStack.cpp +++ b/src/presolve/HighsPostsolveStack.cpp @@ -1353,9 +1353,7 @@ void HighsPostsolveStack::DuplicateColumn::transformToPresolvedSpace( void HighsPostsolveStack::SlackColSubstitution::undo( const HighsOptions& options, const std::vector& rowValues, - const std::vector& colValues, HighsSolution& solution, - HighsBasis& basis) { - + HighsSolution& solution, HighsBasis& basis) { // Taken from HighsPostsolveStack::FreeColSubstitution::undo( // // a (removed) cut may have been used in this reduction. @@ -1378,33 +1376,47 @@ void HighsPostsolveStack::SlackColSubstitution::undo( assert(colCoef != 0); // Row values aren't fully postsolved, so why do this? - if (isModelRow) solution.row_value[row] = + if (isModelRow) + solution.row_value[row] = double(rowValue + colCoef * solution.col_value[col]); - printf("HighsPostsolveStack::SlackColSubstitution::undo rowValue = %g\n", double(rowValue)); + solution.col_value[col] = double((rhs - rowValue) / colCoef); + printf( + "\nHighsPostsolveStack::SlackColSubstitution::undo rowValue = %g; " + "colValue = %g\n", + double(rowValue), solution.col_value[col]); // if no dual values requested, return here if (!solution.dual_valid) return; // compute the row dual value such that reduced cost of basic column is 0 + double save_row_dual = solution.row_dual[row]; if (isModelRow) { solution.row_dual[row] = 0; - HighsCDouble dualval = colCost; - for (const auto& colVal : colValues) { - if (static_cast(colVal.index) < solution.row_dual.size()) - dualval -= colVal.value * solution.row_dual[colVal.index]; - } + HighsCDouble dualval = HighsCDouble(colCost); + dualval = -colCoef * solution.row_dual[row]; solution.row_dual[row] = double(dualval / colCoef); } solution.col_dual[col] = 0; + printf( + "HighsPostsolveStack::SlackColSubstitution::undo OgRowDual = %g; rowDual " + "= %g; colDual = %g\n", + save_row_dual, solution.row_dual[row], solution.col_dual[col]); // set basis status if necessary if (!basis.valid) return; basis.col_status[col] = HighsBasisStatus::kBasic; + HighsBasisStatus save_row_basis_status = basis.row_status[row]; if (isModelRow) - basis.row_status[row] = computeRowStatus(solution.row_dual[row], RowType::kEq); + basis.row_status[row] = + computeRowStatus(solution.row_dual[row], RowType::kEq); + printf( + "HighsPostsolveStack::SlackColSubstitution::undo OgRowStatus = %d; " + "RowStatus = %d; ColStatus = %d\n", + int(save_row_basis_status), int(basis.row_status[row]), + int(basis.col_status[col])); } } // namespace presolve diff --git a/src/presolve/HighsPostsolveStack.h b/src/presolve/HighsPostsolveStack.h index 954ffd2732..418462103c 100644 --- a/src/presolve/HighsPostsolveStack.h +++ b/src/presolve/HighsPostsolveStack.h @@ -222,15 +222,11 @@ class HighsPostsolveStack { struct SlackColSubstitution { double rhs; double colCost; - double colLower; - double colUpper; - // double colCoeff; HighsInt row; HighsInt col; void undo(const HighsOptions& options, - const std::vector& rowValues, - const std::vector& colValues, HighsSolution& solution, + const std::vector& rowValues, HighsSolution& solution, HighsBasis& basis); }; @@ -339,24 +335,17 @@ class HighsPostsolveStack { reductionAdded(ReductionType::kFreeColSubstitution); } - template + template void slackColSubstitution(HighsInt row, HighsInt col, double rhs, - double colCost, double colLower, double colUpper, //double colCoeff, - const HighsMatrixSlice& rowVec, - const HighsMatrixSlice& colVec) { + double colCost, + const HighsMatrixSlice& rowVec) { rowValues.clear(); for (const HighsSliceNonzero& rowVal : rowVec) rowValues.emplace_back(origColIndex[rowVal.index()], rowVal.value()); - colValues.clear(); - for (const HighsSliceNonzero& colVal : colVec) - colValues.emplace_back(origRowIndex[colVal.index()], colVal.value()); - - reductionValues.push(SlackColSubstitution{rhs, colCost, colLower, colUpper, //colCoeff, - origRowIndex[row], - origColIndex[col]}); + reductionValues.push(SlackColSubstitution{rhs, colCost, origRowIndex[row], + origColIndex[col]}); reductionValues.push(rowValues); - reductionValues.push(colValues); reductionAdded(ReductionType::kSlackColSubstitution); } @@ -749,10 +738,9 @@ class HighsPostsolveStack { } case ReductionType::kSlackColSubstitution: { SlackColSubstitution reduction; - reductionValues.pop(colValues); reductionValues.pop(rowValues); reductionValues.pop(reduction); - reduction.undo(options, rowValues, colValues, solution, basis); + reduction.undo(options, rowValues, solution, basis); break; } default: @@ -934,10 +922,9 @@ class HighsPostsolveStack { } case ReductionType::kSlackColSubstitution: { SlackColSubstitution reduction; - reductionValues.pop(colValues); reductionValues.pop(rowValues); reductionValues.pop(reduction); - reduction.undo(options, rowValues, colValues, solution, basis); + reduction.undo(options, rowValues, solution, basis); break; } } diff --git a/src/qpsolver/dantzigpricing.hpp b/src/qpsolver/dantzigpricing.hpp index 5a62dce7fb..3c305529b2 100644 --- a/src/qpsolver/dantzigpricing.hpp +++ b/src/qpsolver/dantzigpricing.hpp @@ -52,7 +52,7 @@ class DantzigPricing : public Pricing { public: DantzigPricing(Runtime& rt, Basis& bas, ReducedCosts& rc) - // clang-format off + // clang-format off : runtime(rt), basis(bas), redcosts(rc) {}; // clang-format on HighsInt price(const QpVector& x, const QpVector& gradient) { diff --git a/src/qpsolver/devexpricing.hpp b/src/qpsolver/devexpricing.hpp index 014d8ce283..6664a98611 100644 --- a/src/qpsolver/devexpricing.hpp +++ b/src/qpsolver/devexpricing.hpp @@ -58,7 +58,7 @@ class DevexPricing : public Pricing { : runtime(rt), basis(bas), redcosts(rc), - // clang-format off + // clang-format off weights(std::vector(rt.instance.num_var, 1.0)) {}; // clang-format on From d0744ed76d717ebc886bf8511fbef477dec1a167 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Fri, 8 Nov 2024 22:40:53 +0000 Subject: [PATCH 7/8] Gets e226 right --- src/presolve/HPresolve.cpp | 2 + src/presolve/HighsPostsolveStack.cpp | 67 +++++++++++++++++----------- src/presolve/HighsPostsolveStack.h | 7 ++- 3 files changed, 48 insertions(+), 28 deletions(-) diff --git a/src/presolve/HPresolve.cpp b/src/presolve/HPresolve.cpp index de5a26204e..8e99f91023 100644 --- a/src/presolve/HPresolve.cpp +++ b/src/presolve/HPresolve.cpp @@ -4429,6 +4429,8 @@ HPresolve::Result HPresolve::removeSlacks( } // postsolve_stack.slackColSubstitution(iRow, iCol, rhs, cost, + lower, + upper, getRowVector(iRow)); markColDeleted(iCol); diff --git a/src/presolve/HighsPostsolveStack.cpp b/src/presolve/HighsPostsolveStack.cpp index 43485006ea..94f7cd0e0c 100644 --- a/src/presolve/HighsPostsolveStack.cpp +++ b/src/presolve/HighsPostsolveStack.cpp @@ -13,6 +13,7 @@ #include #include "lp_data/HConst.h" +#include "lp_data/HighsModelUtils.h" // For debugging #include "lp_data/HighsOptions.h" #include "util/HighsCDouble.h" #include "util/HighsUtils.h" @@ -1381,42 +1382,54 @@ void HighsPostsolveStack::SlackColSubstitution::undo( double(rowValue + colCoef * solution.col_value[col]); solution.col_value[col] = double((rhs - rowValue) / colCoef); + double rowLower = colCoef > 0 ? rhs - colCoef * slackUpper : rhs - colCoef * slackLower; + double rowUpper = colCoef > 0 ? rhs - colCoef * slackLower : rhs - colCoef * slackUpper; + printf( - "\nHighsPostsolveStack::SlackColSubstitution::undo rowValue = %g; " - "colValue = %g\n", - double(rowValue), solution.col_value[col]); + "\nHighsPostsolveStack::SlackColSubstitution::undo rowValue = %11.5g; " + "bounds [%11.5g, %11.5g] colCoef = %11.6g\n", + double(rowValue), rowLower, rowUpper, colCoef); + printf( + "HighsPostsolveStack::SlackColSubstitution::undo colValue = %11.5g, bounds [%11.5g, %11.5g]\n", + solution.col_value[col], slackLower, slackUpper); - // if no dual values requested, return here + // If no dual values requested, return here if (!solution.dual_valid) return; - // compute the row dual value such that reduced cost of basic column is 0 - double save_row_dual = solution.row_dual[row]; + // Row retains its dual value, and column has this dual value scaled by coeff if (isModelRow) { - solution.row_dual[row] = 0; - HighsCDouble dualval = HighsCDouble(colCost); - dualval = -colCoef * solution.row_dual[row]; - solution.row_dual[row] = double(dualval / colCoef); + solution.col_dual[col] = - solution.row_dual[row] / colCoef; + printf( + "HighsPostsolveStack::SlackColSubstitution::undo rowDual = %11.5g; colDual = %11.5g\n", + solution.row_dual[row], solution.col_dual[col]); } - solution.col_dual[col] = 0; - printf( - "HighsPostsolveStack::SlackColSubstitution::undo OgRowDual = %g; rowDual " - "= %g; colDual = %g\n", - save_row_dual, solution.row_dual[row], solution.col_dual[col]); - - // set basis status if necessary + // Set basis status if necessary if (!basis.valid) return; - basis.col_status[col] = HighsBasisStatus::kBasic; - HighsBasisStatus save_row_basis_status = basis.row_status[row]; - if (isModelRow) - basis.row_status[row] = - computeRowStatus(solution.row_dual[row], RowType::kEq); - printf( - "HighsPostsolveStack::SlackColSubstitution::undo OgRowStatus = %d; " - "RowStatus = %d; ColStatus = %d\n", - int(save_row_basis_status), int(basis.row_status[row]), - int(basis.col_status[col])); + // If row is basic, then slack is basic, otherwise row retains its status + if (isModelRow) { + HighsBasisStatus save_row_basis_status = basis.row_status[row]; + if (basis.row_status[row] == HighsBasisStatus::kBasic) { + basis.col_status[col] = HighsBasisStatus::kBasic; + basis.row_status[row] = computeRowStatus(solution.row_dual[row], RowType::kEq); + } else if (basis.row_status[row] == HighsBasisStatus::kLower) { + basis.col_status[col] = colCoef > 0 ? HighsBasisStatus::kUpper : HighsBasisStatus::kLower; + } else { + basis.col_status[col] = colCoef > 0 ? HighsBasisStatus::kLower : HighsBasisStatus::kUpper; + } + printf("HighsPostsolveStack::SlackColSubstitution::undo OgRowStatus = %s; " + "RowStatus = %s; ColStatus = %s\n", + utilBasisStatusToString(save_row_basis_status).c_str(), utilBasisStatusToString(basis.row_status[row]).c_str(), + utilBasisStatusToString(basis.col_status[col]).c_str()); + if (basis.col_status[col] == HighsBasisStatus::kLower) { + assert(solution.col_dual[col] > 0); + } else if (basis.col_status[col] == HighsBasisStatus::kUpper) { + assert(solution.col_dual[col] < 0); + } + } else { + basis.col_status[col] = HighsBasisStatus::kNonbasic; + } } } // namespace presolve diff --git a/src/presolve/HighsPostsolveStack.h b/src/presolve/HighsPostsolveStack.h index 418462103c..8b3f97200e 100644 --- a/src/presolve/HighsPostsolveStack.h +++ b/src/presolve/HighsPostsolveStack.h @@ -222,6 +222,8 @@ class HighsPostsolveStack { struct SlackColSubstitution { double rhs; double colCost; + double slackLower; + double slackUpper; HighsInt row; HighsInt col; @@ -338,12 +340,15 @@ class HighsPostsolveStack { template void slackColSubstitution(HighsInt row, HighsInt col, double rhs, double colCost, + double slackLower, double slackUpper, const HighsMatrixSlice& rowVec) { rowValues.clear(); for (const HighsSliceNonzero& rowVal : rowVec) rowValues.emplace_back(origColIndex[rowVal.index()], rowVal.value()); - reductionValues.push(SlackColSubstitution{rhs, colCost, origRowIndex[row], + reductionValues.push(SlackColSubstitution{rhs, colCost, + slackLower, slackUpper, + origRowIndex[row], origColIndex[col]}); reductionValues.push(rowValues); reductionAdded(ReductionType::kSlackColSubstitution); From f769209b9df6eeb3f86e6464b4778288789c4acd Mon Sep 17 00:00:00 2001 From: JAJHall Date: Sun, 10 Nov 2024 10:31:36 +0000 Subject: [PATCH 8/8] Cleaned up HPresolve::removeSlacks but left off by default since not advantageous on average --- check/TestPresolve.cpp | 2 +- src/lp_data/Highs.cpp | 44 +++--------------------- src/presolve/HPresolve.cpp | 17 ++-------- src/presolve/HighsPostsolveStack.cpp | 50 +++++++++++----------------- src/presolve/HighsPostsolveStack.h | 11 ++---- 5 files changed, 29 insertions(+), 95 deletions(-) diff --git a/check/TestPresolve.cpp b/check/TestPresolve.cpp index 7d2e5d9604..e4263c010e 100644 --- a/check/TestPresolve.cpp +++ b/check/TestPresolve.cpp @@ -616,7 +616,7 @@ TEST_CASE("presolve-slacks", "[highs_test_presolve]") { lp.a_matrix_.index_ = {0, 0}; lp.a_matrix_.value_ = {1, 1}; Highs h; - h.setOptionValue("output_flag", dev_run); + // h.setOptionValue("output_flag", dev_run); REQUIRE(h.passModel(lp) == HighsStatus::kOk); REQUIRE(h.presolve() == HighsStatus::kOk); REQUIRE(h.getPresolvedLp().num_col_ == 0); diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index 894b703aac..f408807367 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -1252,46 +1252,6 @@ HighsStatus Highs::run() { bool have_optimal_solution = false; // ToDo Put solution of presolved problem in a separate method - // if (!this->options_.presolve_remove_slacks) { - HighsLp& reduced_lp = presolve_.getReducedProblem(); - HighsInt num_double_slack = 0; - HighsInt num_slack = 0; - HighsInt num_zero_cost_slack = 0; - HighsInt num_unit_coeff_slack = 0; - double min_slack_coeff = kHighsInf; - double max_slack_coeff = -kHighsInf; - std::vector found_slack; - found_slack.assign(reduced_lp.num_row_, false); - for (HighsInt iCol = 0; iCol < reduced_lp.num_col_; iCol++) { - HighsInt nnz = reduced_lp.a_matrix_.start_[iCol + 1] - - reduced_lp.a_matrix_.start_[iCol]; - if (nnz != 1) continue; - HighsInt iRow = - reduced_lp.a_matrix_.index_[reduced_lp.a_matrix_.start_[iCol]]; - if (found_slack[iRow]) { - num_double_slack++; - continue; - } - if (reduced_lp.row_lower_[iRow] != reduced_lp.row_upper_[iRow]) continue; - num_slack++; - printf("Column %d is slack\n", int(iCol)); - double coeff = std::fabs( - reduced_lp.a_matrix_.value_[reduced_lp.a_matrix_.start_[iCol]]); - if (coeff == 1.0) num_unit_coeff_slack++; - min_slack_coeff = std::min(coeff, min_slack_coeff); - max_slack_coeff = std::max(coeff, max_slack_coeff); - found_slack[iRow] = true; - if (reduced_lp.col_cost_[iCol] == 0) num_zero_cost_slack++; - } - printf( - "grepSlack,model,col,slack,unit coeff,zero_cost,double,min coeff, " - "max_coeff\n"); - printf("grepSlack,%s,%d, %d, %d, %d, %d, %g, %g\n", - this->model_.lp_.model_name_.c_str(), int(reduced_lp.num_col_), - int(num_slack), int(num_unit_coeff_slack), int(num_zero_cost_slack), - int(num_double_slack), min_slack_coeff, max_slack_coeff); - // } - switch (model_presolve_status_) { case HighsPresolveStatus::kNotPresolved: { ekk_instance_.lp_name_ = "Original LP"; @@ -1565,6 +1525,10 @@ HighsStatus Highs::run() { options_ = save_options; if (return_status == HighsStatus::kError) return returnFromRun(return_status, undo_mods); + if (postsolve_iteration_count > 0) + highsLogUser(options_.log_options, HighsLogType::kInfo, + "Required %d simplex iterations after postsolve\n", + int(postsolve_iteration_count)); } } else { highsLogUser(log_options, HighsLogType::kError, diff --git a/src/presolve/HPresolve.cpp b/src/presolve/HPresolve.cpp index 8e99f91023..61d18979f0 100644 --- a/src/presolve/HPresolve.cpp +++ b/src/presolve/HPresolve.cpp @@ -4361,9 +4361,8 @@ HPresolve::Result HPresolve::presolve(HighsPostsolveStack& postsolve_stack) { } // Now consider removing slacks - if (options->presolve_remove_slacks) { + if (options->presolve_remove_slacks) HPRESOLVE_CHECKED_CALL(removeSlacks(postsolve_stack)); - } report(); } else { @@ -4382,10 +4381,8 @@ HPresolve::Result HPresolve::presolve(HighsPostsolveStack& postsolve_stack) { HPresolve::Result HPresolve::removeSlacks( HighsPostsolveStack& postsolve_stack) { - // singletonColumns data structure appears not to be retained + // SingletonColumns data structure appears not to be retained // throughout presolve - // - bool unit_coeff_only = true; for (HighsInt iCol = 0; iCol != model->num_col_; ++iCol) { if (colDeleted[iCol]) continue; if (colsize[iCol] != 1) continue; @@ -4400,11 +4397,6 @@ HPresolve::Result HPresolve::removeSlacks( double cost = model->col_cost_[iCol]; double rhs = model->row_lower_[iRow]; double coeff = Avalue[coliter]; - printf( - "Col %d is continuous and is singleton in equality row %d with cost " - "%g, bounds [%g, %g], coeff %g and RHS = %g\n", - int(iCol), int(iRow), cost, lower, upper, coeff, rhs); - if (unit_coeff_only && std::fabs(coeff) != 1.0) continue; assert(coeff); // Slack is s = (rhs - a^Tx)/coeff // @@ -4428,10 +4420,7 @@ HPresolve::Result HPresolve::removeSlacks( model->offset_ += multiplier * rhs; } // - postsolve_stack.slackColSubstitution(iRow, iCol, rhs, cost, - lower, - upper, - getRowVector(iRow)); + postsolve_stack.slackColSubstitution(iRow, iCol, rhs, getRowVector(iRow)); markColDeleted(iCol); diff --git a/src/presolve/HighsPostsolveStack.cpp b/src/presolve/HighsPostsolveStack.cpp index 94f7cd0e0c..2e78d9a1a7 100644 --- a/src/presolve/HighsPostsolveStack.cpp +++ b/src/presolve/HighsPostsolveStack.cpp @@ -13,7 +13,7 @@ #include #include "lp_data/HConst.h" -#include "lp_data/HighsModelUtils.h" // For debugging +#include "lp_data/HighsModelUtils.h" // For debugging #2001 #include "lp_data/HighsOptions.h" #include "util/HighsCDouble.h" #include "util/HighsUtils.h" @@ -1355,10 +1355,7 @@ void HighsPostsolveStack::DuplicateColumn::transformToPresolvedSpace( void HighsPostsolveStack::SlackColSubstitution::undo( const HighsOptions& options, const std::vector& rowValues, HighsSolution& solution, HighsBasis& basis) { - // Taken from HighsPostsolveStack::FreeColSubstitution::undo( - // - // a (removed) cut may have been used in this reduction. - // + bool debug_print = false; // May have to determine row dual and basis status unless doing // primal-only transformation in MIP solver, in which case row may // no longer exist if it corresponds to a removed cut, so have to @@ -1382,50 +1379,41 @@ void HighsPostsolveStack::SlackColSubstitution::undo( double(rowValue + colCoef * solution.col_value[col]); solution.col_value[col] = double((rhs - rowValue) / colCoef); - double rowLower = colCoef > 0 ? rhs - colCoef * slackUpper : rhs - colCoef * slackLower; - double rowUpper = colCoef > 0 ? rhs - colCoef * slackLower : rhs - colCoef * slackUpper; - - printf( - "\nHighsPostsolveStack::SlackColSubstitution::undo rowValue = %11.5g; " - "bounds [%11.5g, %11.5g] colCoef = %11.6g\n", - double(rowValue), rowLower, rowUpper, colCoef); - printf( - "HighsPostsolveStack::SlackColSubstitution::undo colValue = %11.5g, bounds [%11.5g, %11.5g]\n", - solution.col_value[col], slackLower, slackUpper); // If no dual values requested, return here if (!solution.dual_valid) return; // Row retains its dual value, and column has this dual value scaled by coeff - if (isModelRow) { - solution.col_dual[col] = - solution.row_dual[row] / colCoef; - printf( - "HighsPostsolveStack::SlackColSubstitution::undo rowDual = %11.5g; colDual = %11.5g\n", - solution.row_dual[row], solution.col_dual[col]); - } + if (isModelRow) solution.col_dual[col] = -solution.row_dual[row] / colCoef; // Set basis status if necessary if (!basis.valid) return; - // If row is basic, then slack is basic, otherwise row retains its status + // If row is basic, then slack is basic, otherwise row retains its status if (isModelRow) { HighsBasisStatus save_row_basis_status = basis.row_status[row]; if (basis.row_status[row] == HighsBasisStatus::kBasic) { basis.col_status[col] = HighsBasisStatus::kBasic; - basis.row_status[row] = computeRowStatus(solution.row_dual[row], RowType::kEq); + basis.row_status[row] = + computeRowStatus(solution.row_dual[row], RowType::kEq); } else if (basis.row_status[row] == HighsBasisStatus::kLower) { - basis.col_status[col] = colCoef > 0 ? HighsBasisStatus::kUpper : HighsBasisStatus::kLower; + basis.col_status[col] = + colCoef > 0 ? HighsBasisStatus::kUpper : HighsBasisStatus::kLower; } else { - basis.col_status[col] = colCoef > 0 ? HighsBasisStatus::kLower : HighsBasisStatus::kUpper; + basis.col_status[col] = + colCoef > 0 ? HighsBasisStatus::kLower : HighsBasisStatus::kUpper; } - printf("HighsPostsolveStack::SlackColSubstitution::undo OgRowStatus = %s; " - "RowStatus = %s; ColStatus = %s\n", - utilBasisStatusToString(save_row_basis_status).c_str(), utilBasisStatusToString(basis.row_status[row]).c_str(), - utilBasisStatusToString(basis.col_status[col]).c_str()); + if (debug_print) + printf( + "HighsPostsolveStack::SlackColSubstitution::undo OgRowStatus = %s; " + "RowStatus = %s; ColStatus = %s\n", + utilBasisStatusToString(save_row_basis_status).c_str(), + utilBasisStatusToString(basis.row_status[row]).c_str(), + utilBasisStatusToString(basis.col_status[col]).c_str()); if (basis.col_status[col] == HighsBasisStatus::kLower) { - assert(solution.col_dual[col] > 0); + assert(solution.col_dual[col] > -options.dual_feasibility_tolerance); } else if (basis.col_status[col] == HighsBasisStatus::kUpper) { - assert(solution.col_dual[col] < 0); + assert(solution.col_dual[col] < options.dual_feasibility_tolerance); } } else { basis.col_status[col] = HighsBasisStatus::kNonbasic; diff --git a/src/presolve/HighsPostsolveStack.h b/src/presolve/HighsPostsolveStack.h index 8b3f97200e..f3409d074c 100644 --- a/src/presolve/HighsPostsolveStack.h +++ b/src/presolve/HighsPostsolveStack.h @@ -221,9 +221,6 @@ class HighsPostsolveStack { struct SlackColSubstitution { double rhs; - double colCost; - double slackLower; - double slackUpper; HighsInt row; HighsInt col; @@ -339,17 +336,13 @@ class HighsPostsolveStack { template void slackColSubstitution(HighsInt row, HighsInt col, double rhs, - double colCost, - double slackLower, double slackUpper, const HighsMatrixSlice& rowVec) { rowValues.clear(); for (const HighsSliceNonzero& rowVal : rowVec) rowValues.emplace_back(origColIndex[rowVal.index()], rowVal.value()); - reductionValues.push(SlackColSubstitution{rhs, colCost, - slackLower, slackUpper, - origRowIndex[row], - origColIndex[col]}); + reductionValues.push( + SlackColSubstitution{rhs, origRowIndex[row], origColIndex[col]}); reductionValues.push(rowValues); reductionAdded(ReductionType::kSlackColSubstitution); }