From 6f9dc9205e6a23e4430d74cbf7d049a05e2ef4e9 Mon Sep 17 00:00:00 2001 From: fwesselm Date: Fri, 10 Jan 2025 11:00:54 +0100 Subject: [PATCH 01/11] Remove usage of std::vector::at from HEkkPrimal.cpp --- src/simplex/HEkkPrimal.cpp | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/simplex/HEkkPrimal.cpp b/src/simplex/HEkkPrimal.cpp index 1c39033e7e..e077da8040 100644 --- a/src/simplex/HEkkPrimal.cpp +++ b/src/simplex/HEkkPrimal.cpp @@ -1212,11 +1212,11 @@ void HEkkPrimal::phase1ChooseRow() { analysis->simplexTimerStart(Chuzr2Clock); pdqsort(ph1SorterR.begin(), ph1SorterR.end()); - double dMaxTheta = ph1SorterR.at(0).first; + double dMaxTheta = ph1SorterR[0].first; double dGradient = fabs(theta_dual); for (size_t i = 0; i < ph1SorterR.size(); i++) { - double dMyTheta = ph1SorterR.at(i).first; - HighsInt index = ph1SorterR.at(i).second; + double dMyTheta = ph1SorterR[i].first; + HighsInt index = ph1SorterR[i].second; HighsInt iRow = index >= 0 ? index : index + num_row; dGradient -= fabs(col_aq.array[iRow]); // Stop when the gradient start to decrease @@ -1231,8 +1231,8 @@ void HEkkPrimal::phase1ChooseRow() { double dMaxAlpha = 0.0; size_t iLast = ph1SorterT.size(); for (size_t i = 0; i < ph1SorterT.size(); i++) { - double dMyTheta = ph1SorterT.at(i).first; - HighsInt index = ph1SorterT.at(i).second; + double dMyTheta = ph1SorterT[i].first; + HighsInt index = ph1SorterT[i].second; HighsInt iRow = index >= 0 ? index : index + num_row; double dAbsAlpha = fabs(col_aq.array[iRow]); // Stop when the theta is too large @@ -1251,7 +1251,7 @@ void HEkkPrimal::phase1ChooseRow() { variable_out = -1; move_out = 0; for (size_t i = iLast; i > 0; i--) { - HighsInt index = ph1SorterT.at(i - 1).second; + HighsInt index = ph1SorterT[i - 1].second; HighsInt iRow = index >= 0 ? index : index + num_row; double dAbsAlpha = fabs(col_aq.array[iRow]); if (dAbsAlpha > dMaxAlpha * 0.1) { From 9e19d1fb02f22f59869769110050043b121a32ec Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 13:21:05 +0000 Subject: [PATCH 02/11] Now bailing out of deleteColsInterface and deleteRowsInterface if no cols/rows are deleted --- check/TestLpModification.cpp | 39 ++++++++++++++++++++++++++++++++++ src/lp_data/HStruct.h | 2 ++ src/lp_data/HighsInterface.cpp | 34 +++++++++++++++++------------ 3 files changed, 61 insertions(+), 14 deletions(-) diff --git a/check/TestLpModification.cpp b/check/TestLpModification.cpp index 30bdac0f06..5c58faa7c7 100644 --- a/check/TestLpModification.cpp +++ b/check/TestLpModification.cpp @@ -2273,4 +2273,43 @@ TEST_CASE("row-wise-get-row-avgas", "[highs_data]") { h.ensureColwise(); testAvgasGetRow(h); testAvgasGetCol(h); + +} + +TEST_CASE("hot-start-after-delete", "[highs_data]") { + Highs h; + h.setOptionValue("output_flag", dev_run); + const HighsInfo& info = h.getInfo(); + const HighsBasis& basis = h.getBasis(); + const HighsSolution& solution = h.getSolution(); + std::string model = "avgas"; + std::string model_file = + std::string(HIGHS_DIR) + "/check/instances/" + model + ".mps"; + h.readModel(model_file); + h.run(); + printf("Initial solve takes %d iterations and yields objective = %g\n", + int(info.simplex_iteration_count), info.objective_function_value); + h.writeSolution("", 1); + double cost; + double lower; + double upper; + HighsInt nnz; + std::vector start(1); + std::vector index(h.getNumRow()); + std::vector value(h.getNumRow()); + HighsInt get_num; + + HighsInt use_col = 1; + printf("Deleting and adding column %1d with status \"%s\" and value %g\n", + int(use_col), h.basisStatusToString(basis.col_status[use_col]).c_str(), solution.col_value[use_col]); + + h.getCols(use_col, use_col, get_num, &cost, &lower, &upper, nnz, start.data(), index.data(), value.data()); + + h.deleteCols(use_col, use_col); + + h.addCol(cost, lower, upper, nnz, index.data(), value.data()); + + h.run(); + printf("After deleting and adding column %1d, solve takes %d iterations and yields objective = %g\n", + int(use_col), int(info.simplex_iteration_count), info.objective_function_value); } diff --git a/src/lp_data/HStruct.h b/src/lp_data/HStruct.h index 77b9e6e11b..020c3bda1d 100644 --- a/src/lp_data/HStruct.h +++ b/src/lp_data/HStruct.h @@ -66,6 +66,8 @@ struct HighsBasis { std::string debug_origin_name = "None"; std::vector col_status; std::vector row_status; + void deleteCols(const HighsIndexCollection& index_collection); + void deleteRows(const HighsIndexCollection& index_collection); void invalidate(); void clear(); }; diff --git a/src/lp_data/HighsInterface.cpp b/src/lp_data/HighsInterface.cpp index e522e46db9..a474d7c8d9 100644 --- a/src/lp_data/HighsInterface.cpp +++ b/src/lp_data/HighsInterface.cpp @@ -587,13 +587,16 @@ void Highs::deleteColsInterface(HighsIndexCollection& index_collection) { lp.deleteCols(index_collection); model_.hessian_.deleteCols(index_collection); - assert(lp.num_col_ <= original_num_col); - if (lp.num_col_ < original_num_col) { - // Nontrivial deletion so reset the model_status and invalidate - // the Highs basis - model_status_ = HighsModelStatus::kNotset; - basis.valid = false; - } + // Bail out if no columns were actually deleted + if (lp.num_col_ == original_num_col) return; + + assert(lp.num_col_ < original_num_col); + + // Nontrivial deletion so reset the model_status and invalidate + // the Highs basis + model_status_ = HighsModelStatus::kNotset; + basis.valid = false; + if (lp.scale_.has_scaling) { deleteScale(lp.scale_.col, index_collection); lp.scale_.col.resize(lp.num_col_); @@ -632,13 +635,16 @@ void Highs::deleteRowsInterface(HighsIndexCollection& index_collection) { HighsInt original_num_row = lp.num_row_; lp.deleteRows(index_collection); - assert(lp.num_row_ <= original_num_row); - if (lp.num_row_ < original_num_row) { - // Nontrivial deletion so reset the model_status and invalidate - // the Highs basis - model_status_ = HighsModelStatus::kNotset; - basis.valid = false; - } + // Bail out if no rows were actually deleted + if (lp.num_row_ == original_num_row) return; + + assert(lp.num_row_ < original_num_row); + + // Nontrivial deletion so reset the model_status and invalidate + // the Highs basis + model_status_ = HighsModelStatus::kNotset; + basis.valid = false; + if (lp.scale_.has_scaling) { deleteScale(lp.scale_.row, index_collection); lp.scale_.row.resize(lp.num_row_); From bdabadf261d53e43d83d133389db4280da269638 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 15:50:32 +0000 Subject: [PATCH 03/11] Introduce HighsBasis::useful --- check/TestLpModification.cpp | 29 +++++++- src/lp_data/HStruct.h | 4 +- src/lp_data/Highs.cpp | 9 +++ src/lp_data/HighsInterface.cpp | 121 +++++++++++++++++++++++++++++---- src/lp_data/HighsSolution.cpp | 17 +++++ 5 files changed, 163 insertions(+), 17 deletions(-) diff --git a/check/TestLpModification.cpp b/check/TestLpModification.cpp index 5c58faa7c7..be24bf6605 100644 --- a/check/TestLpModification.cpp +++ b/check/TestLpModification.cpp @@ -2279,6 +2279,7 @@ TEST_CASE("row-wise-get-row-avgas", "[highs_data]") { TEST_CASE("hot-start-after-delete", "[highs_data]") { Highs h; h.setOptionValue("output_flag", dev_run); + const HighsLp& lp = h.getLp(); const HighsInfo& info = h.getInfo(); const HighsBasis& basis = h.getBasis(); const HighsSolution& solution = h.getSolution(); @@ -2290,26 +2291,48 @@ TEST_CASE("hot-start-after-delete", "[highs_data]") { printf("Initial solve takes %d iterations and yields objective = %g\n", int(info.simplex_iteration_count), info.objective_function_value); h.writeSolution("", 1); + for (HighsInt iCol = 0; iCol < lp.num_col_; iCol++) + printf("Col %2d is %s\n", int(iCol), basis.col_status[iCol] == HighsBasisStatus::kBasic ? "basic" : "nonbasic"); + printf("\n"); + for (HighsInt iRow = 0; iRow < lp.num_row_; iRow++) + printf("Row %2d is %s\n", int(iRow), basis.row_status[iRow] == HighsBasisStatus::kBasic ? "basic" : "nonbasic"); + double cost; double lower; double upper; HighsInt nnz; std::vector start(1); - std::vector index(h.getNumRow()); - std::vector value(h.getNumRow()); + std::vector index(lp.num_row_); + std::vector value(lp.num_row_); HighsInt get_num; HighsInt use_col = 1; - printf("Deleting and adding column %1d with status \"%s\" and value %g\n", + printf("\nDeleting and adding column %1d with status \"%s\" and value %g\n", int(use_col), h.basisStatusToString(basis.col_status[use_col]).c_str(), solution.col_value[use_col]); h.getCols(use_col, use_col, get_num, &cost, &lower, &upper, nnz, start.data(), index.data(), value.data()); h.deleteCols(use_col, use_col); + basis.printScalars(); h.addCol(cost, lower, upper, nnz, index.data(), value.data()); h.run(); printf("After deleting and adding column %1d, solve takes %d iterations and yields objective = %g\n", int(use_col), int(info.simplex_iteration_count), info.objective_function_value); + + HighsInt use_row = 1; + printf("\nDeleting and adding row %1d with status \"%s\" and value %g\n", + int(use_row), h.basisStatusToString(basis.row_status[use_row]).c_str(), solution.row_value[use_row]); + + h.getRows(use_row, use_row, get_num, &lower, &upper, nnz, start.data(), index.data(), value.data()); + + h.deleteRows(use_row, use_row); + basis.printScalars(); + + h.addRow(lower, upper, nnz, index.data(), value.data()); + + h.run(); + printf("After deleting and adding row %1d, solve takes %d iterations and yields objective = %g\n", + int(use_row), int(info.simplex_iteration_count), info.objective_function_value); } diff --git a/src/lp_data/HStruct.h b/src/lp_data/HStruct.h index 020c3bda1d..51d55acc46 100644 --- a/src/lp_data/HStruct.h +++ b/src/lp_data/HStruct.h @@ -66,8 +66,8 @@ struct HighsBasis { std::string debug_origin_name = "None"; std::vector col_status; std::vector row_status; - void deleteCols(const HighsIndexCollection& index_collection); - void deleteRows(const HighsIndexCollection& index_collection); + void print() const; + void printScalars() const; void invalidate(); void clear(); }; diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index df38dad281..672c3a2bcf 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -1265,6 +1265,15 @@ HighsStatus Highs::solve() { // is available, simplex should surely be chosen. const bool solver_will_use_basis = options_.solver == kSimplexString || options_.solver == kHighsChooseString; + const bool has_basis = + basis_.col_status.size() == static_cast(model_.lp_.num_col_) && + basis_.row_status.size() == static_cast(model_.lp_.num_row_); + if (has_basis && !basis_.valid) { + basis_.print(); + printf("Highs::solve() has_basis && !basis_.valid\n"); + assert(111==123); + } + if ((basis_.valid || options_.presolve == kHighsOffString || unconstrained_lp) && solver_will_use_basis) { diff --git a/src/lp_data/HighsInterface.cpp b/src/lp_data/HighsInterface.cpp index a474d7c8d9..743f7dd3f7 100644 --- a/src/lp_data/HighsInterface.cpp +++ b/src/lp_data/HighsInterface.cpp @@ -576,6 +576,76 @@ HighsStatus Highs::addRowsInterface(HighsInt ext_num_new_row, return return_status; } +void deleteBasisEntries(std::vector& status, + bool& deleted_basic, + bool& deleted_nonbasic, + const HighsIndexCollection& index_collection, + const HighsInt entry_dim) { + assert(ok(index_collection)); + assert(static_cast(entry_dim) == status.size()); + HighsInt from_k; + HighsInt to_k; + limits(index_collection, from_k, to_k); + if (from_k > to_k) return; + + HighsInt delete_from_entry; + HighsInt delete_to_entry; + HighsInt keep_from_entry; + HighsInt keep_to_entry = -1; + HighsInt current_set_entry = 0; + HighsInt new_num_entry = 0; + deleted_basic = false; + deleted_nonbasic = false; + for (HighsInt k = from_k; k <= to_k; k++) { + updateOutInIndex(index_collection, delete_from_entry, delete_to_entry, + keep_from_entry, keep_to_entry, current_set_entry); + // Account for the initial entries being kept + if (k == from_k) new_num_entry = delete_from_entry; + if (delete_to_entry >= entry_dim - 1) break; + assert(delete_to_entry < entry_dim); + // Identify whether a basic or a nonbasic entry has been deleted + for (HighsInt entry = delete_from_entry; entry <= delete_to_entry; entry++) { + if (status[entry] == HighsBasisStatus::kBasic) { + deleted_basic = true; + } else { + deleted_nonbasic = true; + } + } + for (HighsInt entry = keep_from_entry; entry <= keep_to_entry; entry++) { + status[new_num_entry] = status[entry]; + new_num_entry++; + } + if (keep_to_entry >= entry_dim - 1) break; + } + status.resize(new_num_entry); +} + +void deleteBasisCols(HighsBasis& basis, + const HighsIndexCollection& index_collection, + const HighsInt original_num_col) { + bool deleted_basic; + bool deleted_nonbasic; + deleteBasisEntries(basis.col_status, + deleted_basic, + deleted_nonbasic, + index_collection, original_num_col); + if (deleted_basic) + basis.valid = false; +} + +void deleteBasisRows(HighsBasis& basis, + const HighsIndexCollection& index_collection, + const HighsInt original_num_row) { + bool deleted_basic; + bool deleted_nonbasic; + deleteBasisEntries(basis.row_status, + deleted_basic, + deleted_nonbasic, + index_collection, original_num_row); + if (deleted_nonbasic) + basis.valid = false; +} + void Highs::deleteColsInterface(HighsIndexCollection& index_collection) { HighsLp& lp = model_.lp_; HighsBasis& basis = basis_; @@ -592,11 +662,19 @@ void Highs::deleteColsInterface(HighsIndexCollection& index_collection) { assert(lp.num_col_ < original_num_col); - // Nontrivial deletion so reset the model_status and invalidate - // the Highs basis + // Nontrivial deletion so reset the model_status and update any + // Highs basis model_status_ = HighsModelStatus::kNotset; + if (basis_.col_status.size() == static_cast(original_num_col)) { + // Have a full set of column basis status values, so maintain + // them, and only invalidate the basis if a basic column has been + // deleted + deleteBasisCols(basis_, index_collection, original_num_col); + } else { + assert(!basis.valid); + } basis.valid = false; - + if (lp.scale_.has_scaling) { deleteScale(lp.scale_.col, index_collection); lp.scale_.col.resize(lp.num_col_); @@ -640,9 +718,17 @@ void Highs::deleteRowsInterface(HighsIndexCollection& index_collection) { assert(lp.num_row_ < original_num_row); - // Nontrivial deletion so reset the model_status and invalidate - // the Highs basis + // Nontrivial deletion so reset the model_status and update any + // Highs basis model_status_ = HighsModelStatus::kNotset; + if (basis_.row_status.size() == static_cast(original_num_row)) { + // Have a full set of row basis status values, so maintain them, + // and only invalidate the basis if a nonbasic row has been + // deleted + deleteBasisRows(basis_, index_collection, original_num_row); + } else { + assert(!basis.valid); + } basis.valid = false; if (lp.scale_.has_scaling) { @@ -1175,17 +1261,21 @@ void Highs::setNonbasicStatusInterface( } void Highs::appendNonbasicColsToBasisInterface(const HighsInt ext_num_new_col) { + if (ext_num_new_col == 0) return; HighsBasis& highs_basis = basis_; if (!highs_basis.valid) return; const bool has_simplex_basis = ekk_instance_.status_.has_basis; SimplexBasis& simplex_basis = ekk_instance_.basis_; HighsLp& lp = model_.lp_; + const bool has_highs_basis = highs_basis.valid && + highs_basis.col_status.size() == static_cast(lp.num_col_) && + highs_basis.row_status.size() == static_cast(lp.num_row_); + if (!has_highs_basis && !has_simplex_basis) return; // Add nonbasic structurals - if (ext_num_new_col == 0) return; HighsInt newNumCol = lp.num_col_ + ext_num_new_col; HighsInt newNumTot = newNumCol + lp.num_row_; - highs_basis.col_status.resize(newNumCol); + if (has_highs_basis) highs_basis.col_status.resize(newNumCol); if (has_simplex_basis) { simplex_basis.nonbasicFlag_.resize(newNumTot); simplex_basis.nonbasicMove_.resize(newNumTot); @@ -1239,7 +1329,7 @@ void Highs::appendNonbasicColsToBasisInterface(const HighsInt ext_num_new_col) { move = kNonbasicMoveZe; } assert(status != HighsBasisStatus::kNonbasic); - highs_basis.col_status[iCol] = status; + if (has_highs_basis) highs_basis.col_status[iCol] = status; if (has_simplex_basis) { assert(move != kIllegalMoveValue); simplex_basis.nonbasicFlag_[iCol] = kNonbasicFlagTrue; @@ -1249,18 +1339,25 @@ void Highs::appendNonbasicColsToBasisInterface(const HighsInt ext_num_new_col) { } void Highs::appendBasicRowsToBasisInterface(const HighsInt ext_num_new_row) { + if (ext_num_new_row == 0) return; HighsBasis& highs_basis = basis_; if (!highs_basis.valid) return; const bool has_simplex_basis = ekk_instance_.status_.has_basis; SimplexBasis& simplex_basis = ekk_instance_.basis_; HighsLp& lp = model_.lp_; + const bool has_highs_basis = highs_basis.valid && + highs_basis.col_status.size() == static_cast(lp.num_col_) && + highs_basis.row_status.size() == static_cast(lp.num_row_); + if (!has_highs_basis && !has_simplex_basis) return; + // Add basic logicals - if (ext_num_new_row == 0) return; // Add the new rows to the Highs basis HighsInt newNumRow = lp.num_row_ + ext_num_new_row; - highs_basis.row_status.resize(newNumRow); - for (HighsInt iRow = lp.num_row_; iRow < newNumRow; iRow++) - highs_basis.row_status[iRow] = HighsBasisStatus::kBasic; + if (has_highs_basis) { + highs_basis.row_status.resize(newNumRow); + for (HighsInt iRow = lp.num_row_; iRow < newNumRow; iRow++) + highs_basis.row_status[iRow] = HighsBasisStatus::kBasic; + } if (has_simplex_basis) { // Add the new rows to the simplex basis HighsInt newNumTot = lp.num_col_ + newNumRow; diff --git a/src/lp_data/HighsSolution.cpp b/src/lp_data/HighsSolution.cpp index 5375a3b6b9..d63aed8964 100644 --- a/src/lp_data/HighsSolution.cpp +++ b/src/lp_data/HighsSolution.cpp @@ -1605,6 +1605,23 @@ void HighsSolution::clear() { void HighsObjectiveSolution::clear() { this->col_value.clear(); } +void HighsBasis::print() const { + this->printScalars(); + for (HighsInt iCol = 0; iCol < HighsInt(this->col_status.size()); iCol++) + printf("Basis: col_status[%2d] = %d\n", int(iCol), int(this->col_status[iCol])); + for (HighsInt iRow = 0; iRow < HighsInt(this->row_status.size()); iRow++) + printf("Basis: row_status[%2d] = %d\n", int(iRow), int(this->row_status[iRow])); +} + +void HighsBasis::printScalars() const { + printf("Basis: valid = %d\n", this->valid); + printf("Basis: alien = %d\n", this->alien); + printf("Basis: was_alien = %d\n", this->was_alien); + printf("Basis: debug_id = %d\n", int(this->debug_id)); + printf("Basis: debug_update_count = %d\n", int(this->debug_update_count)); + printf("Basis: debug_origin_name = %s\n", this->debug_origin_name.c_str()); +} + void HighsBasis::invalidate() { this->valid = false; this->alien = true; From 07faec298521550e57586cf7f48156a09a3c7fd1 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 16:35:39 +0000 Subject: [PATCH 04/11] Introduced HighsBasis::useful and all unt tests pass --- src/lp_data/HStruct.h | 1 + src/lp_data/Highs.cpp | 16 +++++++++++++--- src/lp_data/HighsInterface.cpp | 26 ++++++++++++++++++-------- src/lp_data/HighsSolution.cpp | 4 ++++ src/lp_data/HighsSolve.cpp | 1 + src/mip/HighsMipSolverData.cpp | 3 +++ src/presolve/HPresolve.cpp | 1 + src/qpsolver/a_quass.cpp | 1 + src/simplex/HEkk.cpp | 1 + 9 files changed, 43 insertions(+), 11 deletions(-) diff --git a/src/lp_data/HStruct.h b/src/lp_data/HStruct.h index 51d55acc46..38e78855df 100644 --- a/src/lp_data/HStruct.h +++ b/src/lp_data/HStruct.h @@ -60,6 +60,7 @@ struct HotStart { struct HighsBasis { bool valid = false; bool alien = true; + bool useful = false; bool was_alien = true; HighsInt debug_id = -1; HighsInt debug_update_count = -1; diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index 672c3a2bcf..3e45df3443 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -716,6 +716,7 @@ HighsStatus Highs::readBasis(const std::string& filename) { // Update the HiGHS basis and invalidate any simplex basis for the model basis_ = read_basis; basis_.valid = true; + basis_.useful = true; // Follow implications of a new HiGHS basis newHighsBasis(); // Can't use returnFromHighs since... @@ -1265,9 +1266,12 @@ HighsStatus Highs::solve() { // is available, simplex should surely be chosen. const bool solver_will_use_basis = options_.solver == kSimplexString || options_.solver == kHighsChooseString; - const bool has_basis = - basis_.col_status.size() == static_cast(model_.lp_.num_col_) && - basis_.row_status.size() == static_cast(model_.lp_.num_row_); + const bool has_basis = basis_.valid && basis_.useful; + if (has_basis) { + assert(basis_.col_status.size() == static_cast(incumbent_lp.num_col_)); + assert(basis_.row_status.size() == static_cast(incumbent_lp.num_row_)); + } + if (basis_.valid) assert(basis_.useful); if (has_basis && !basis_.valid) { basis_.print(); printf("Highs::solve() has_basis && !basis_.valid\n"); @@ -1425,6 +1429,7 @@ HighsStatus Highs::solve() { basis_.debug_origin_name = "Presolve to empty"; basis_.valid = true; basis_.alien = false; + basis_.useful = true; basis_.was_alien = false; solution_.value_valid = true; solution_.dual_valid = true; @@ -1560,6 +1565,7 @@ HighsStatus Highs::solve() { solution_.dual_valid = true; // Set basis and its status basis_.valid = true; + basis_.useful = true; basis_.col_status = presolve_.data_.recovered_basis_.col_status; basis_.row_status = presolve_.data_.recovered_basis_.row_status; basis_.debug_origin_name += ": after postsolve"; @@ -2348,6 +2354,7 @@ HighsStatus Highs::setBasis(const HighsBasis& basis, basis_ = basis; } basis_.valid = true; + basis_.useful = true; if (origin != "") basis_.debug_origin_name = origin; assert(basis_.debug_origin_name != ""); assert(!basis_.alien); @@ -4166,6 +4173,7 @@ HighsStatus Highs::callRunPostsolve(const HighsSolution& solution, // Set basis and its status // // basis_.valid = true; + // basis_.useful = true; // basis_.col_status = presolve_.data_.recovered_basis_.col_status; // basis_.row_status = presolve_.data_.recovered_basis_.row_status; basis_.debug_origin_name += ": after postsolve"; @@ -4279,10 +4287,12 @@ void Highs::forceHighsSolutionBasisSize() { if (basis_.col_status.size() != static_cast(model_.lp_.num_col_)) { basis_.col_status.resize(model_.lp_.num_col_); basis_.valid = false; + basis_.useful = false; } if (basis_.row_status.size() != static_cast(model_.lp_.num_row_)) { basis_.row_status.resize(model_.lp_.num_row_); basis_.valid = false; + basis_.useful = false; } } diff --git a/src/lp_data/HighsInterface.cpp b/src/lp_data/HighsInterface.cpp index 743f7dd3f7..31eb3edcd2 100644 --- a/src/lp_data/HighsInterface.cpp +++ b/src/lp_data/HighsInterface.cpp @@ -665,7 +665,8 @@ void Highs::deleteColsInterface(HighsIndexCollection& index_collection) { // Nontrivial deletion so reset the model_status and update any // Highs basis model_status_ = HighsModelStatus::kNotset; - if (basis_.col_status.size() == static_cast(original_num_col)) { + if (basis_.useful) { + assert(basis_.col_status.size() == static_cast(original_num_col)); // Have a full set of column basis status values, so maintain // them, and only invalidate the basis if a basic column has been // deleted @@ -721,7 +722,8 @@ void Highs::deleteRowsInterface(HighsIndexCollection& index_collection) { // Nontrivial deletion so reset the model_status and update any // Highs basis model_status_ = HighsModelStatus::kNotset; - if (basis_.row_status.size() == static_cast(original_num_row)) { + if (basis_.useful) { + assert(basis_.row_status.size() == static_cast(original_num_row)); // Have a full set of row basis status values, so maintain them, // and only invalidate the basis if a nonbasic row has been // deleted @@ -1267,9 +1269,12 @@ void Highs::appendNonbasicColsToBasisInterface(const HighsInt ext_num_new_col) { const bool has_simplex_basis = ekk_instance_.status_.has_basis; SimplexBasis& simplex_basis = ekk_instance_.basis_; HighsLp& lp = model_.lp_; - const bool has_highs_basis = highs_basis.valid && - highs_basis.col_status.size() == static_cast(lp.num_col_) && - highs_basis.row_status.size() == static_cast(lp.num_row_); + const bool has_highs_basis = highs_basis.valid && highs_basis.useful; + + if (has_highs_basis) { + assert(highs_basis.col_status.size() == static_cast(lp.num_col_)); + assert(highs_basis.row_status.size() == static_cast(lp.num_row_)); + } if (!has_highs_basis && !has_simplex_basis) return; // Add nonbasic structurals @@ -1345,9 +1350,13 @@ void Highs::appendBasicRowsToBasisInterface(const HighsInt ext_num_new_row) { const bool has_simplex_basis = ekk_instance_.status_.has_basis; SimplexBasis& simplex_basis = ekk_instance_.basis_; HighsLp& lp = model_.lp_; - const bool has_highs_basis = highs_basis.valid && - highs_basis.col_status.size() == static_cast(lp.num_col_) && - highs_basis.row_status.size() == static_cast(lp.num_row_); + const bool has_highs_basis = highs_basis.valid && highs_basis.useful; + + if (has_highs_basis) { + assert(highs_basis.col_status.size() == static_cast(lp.num_col_)); + assert(highs_basis.row_status.size() == static_cast(lp.num_row_)); + } + if (!has_highs_basis && !has_simplex_basis) return; // Add basic logicals @@ -1641,6 +1650,7 @@ HighsStatus Highs::setHotStartInterface(const HotStart& hot_start) { nonbasicMove[num_col + iRow] = move; } basis_.valid = true; + basis_.useful = true; ekk_instance_.status_.has_basis = true; ekk_instance_.setNlaRefactorInfo(); ekk_instance_.updateStatus(LpAction::kHotStart); diff --git a/src/lp_data/HighsSolution.cpp b/src/lp_data/HighsSolution.cpp index d63aed8964..fb08dc9c2f 100644 --- a/src/lp_data/HighsSolution.cpp +++ b/src/lp_data/HighsSolution.cpp @@ -1319,6 +1319,7 @@ HighsStatus ipxBasicSolutionToHighsBasicSolution( highs_solution.value_valid = true; highs_solution.dual_valid = true; highs_basis.valid = true; + highs_basis.useful = true; return HighsStatus::kOk; } @@ -1606,6 +1607,7 @@ void HighsSolution::clear() { void HighsObjectiveSolution::clear() { this->col_value.clear(); } void HighsBasis::print() const { + if (!this->useful) return; this->printScalars(); for (HighsInt iCol = 0; iCol < HighsInt(this->col_status.size()); iCol++) printf("Basis: col_status[%2d] = %d\n", int(iCol), int(this->col_status[iCol])); @@ -1616,6 +1618,7 @@ void HighsBasis::print() const { void HighsBasis::printScalars() const { printf("Basis: valid = %d\n", this->valid); printf("Basis: alien = %d\n", this->alien); + printf("Basis: useful = %d\n", this->useful); printf("Basis: was_alien = %d\n", this->was_alien); printf("Basis: debug_id = %d\n", int(this->debug_id)); printf("Basis: debug_update_count = %d\n", int(this->debug_update_count)); @@ -1625,6 +1628,7 @@ void HighsBasis::printScalars() const { void HighsBasis::invalidate() { this->valid = false; this->alien = true; + this->useful = false; this->was_alien = true; this->debug_id = -1; this->debug_update_count = -1; diff --git a/src/lp_data/HighsSolve.cpp b/src/lp_data/HighsSolve.cpp index 1062b35e5f..acc81abf43 100644 --- a/src/lp_data/HighsSolve.cpp +++ b/src/lp_data/HighsSolve.cpp @@ -390,6 +390,7 @@ HighsStatus solveUnconstrainedLp(const HighsOptions& options, const HighsLp& lp, solution.value_valid = true; solution.dual_valid = true; basis.valid = true; + basis.useful = true; highs_info.basis_validity = kBasisValidityValid; setSolutionStatus(highs_info); if (highs_info.num_primal_infeasibilities) { diff --git a/src/mip/HighsMipSolverData.cpp b/src/mip/HighsMipSolverData.cpp index 928897c980..4c411f8ffc 100644 --- a/src/mip/HighsMipSolverData.cpp +++ b/src/mip/HighsMipSolverData.cpp @@ -1265,6 +1265,7 @@ void HighsMipSolverData::performRestart() { root_basis.row_status.resize(postSolveStack.getOrigNumRow(), HighsBasisStatus::kBasic); root_basis.valid = true; + root_basis.useful = true; for (HighsInt i = 0; i < mipsolver.model_->num_col_; ++i) root_basis.col_status[postSolveStack.getOrigColIndex(i)] = @@ -1379,6 +1380,7 @@ void HighsMipSolverData::basisTransfer() { firstrootbasis.row_status.assign(numRow, HighsBasisStatus::kNonbasic); firstrootbasis.valid = true; firstrootbasis.alien = true; + firstrootbasis.useful = true; for (HighsInt i = 0; i < numRow; ++i) { HighsBasisStatus status = @@ -1929,6 +1931,7 @@ void HighsMipSolverData::evaluateRootNode() { firstrootbasis.row_status.assign(mipsolver.numRow(), HighsBasisStatus::kBasic); firstrootbasis.valid = true; + firstrootbasis.useful = true; } if (cutpool.getNumCuts() != 0) { diff --git a/src/presolve/HPresolve.cpp b/src/presolve/HPresolve.cpp index eb075e53fb..dbb7854d94 100644 --- a/src/presolve/HPresolve.cpp +++ b/src/presolve/HPresolve.cpp @@ -6330,6 +6330,7 @@ void HPresolve::debug(const HighsLp& lp, const HighsOptions& options) { temp_sol.row_value.resize(model.num_row_); calculateRowValuesQuad(model, sol); temp_basis.valid = true; + temp_basis.useful = true; refineBasis(model, temp_sol, temp_basis); Highs highs; highs.passOptions(options); diff --git a/src/qpsolver/a_quass.cpp b/src/qpsolver/a_quass.cpp index 56ff800596..d60bc999e5 100644 --- a/src/qpsolver/a_quass.cpp +++ b/src/qpsolver/a_quass.cpp @@ -116,6 +116,7 @@ static QpAsmStatus quass2highs(Instance& instance, Settings& settings, } highs_basis.valid = true; highs_basis.alien = false; + highs_basis.useful = true; return qp_asm_return_status; } diff --git a/src/simplex/HEkk.cpp b/src/simplex/HEkk.cpp index 2dd2ff85ae..c6da49da49 100644 --- a/src/simplex/HEkk.cpp +++ b/src/simplex/HEkk.cpp @@ -1465,6 +1465,7 @@ HighsBasis HEkk::getHighsBasis(HighsLp& use_lp) const { } highs_basis.valid = true; highs_basis.alien = false; + highs_basis.useful = true; highs_basis.was_alien = false; highs_basis.debug_id = (HighsInt)(build_synthetic_tick_ + total_synthetic_tick_); From b03608c811f4ea329b60550df86477d2745ca3cf Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 17:50:34 +0000 Subject: [PATCH 05/11] All unit tests pass --- src/lp_data/HStruct.h | 14 ++++++++++++-- src/lp_data/Highs.cpp | 9 ++------- src/lp_data/HighsSolution.cpp | 33 ++++++++++++++++++++------------- src/simplex/HApp.h | 12 ++++++++++++ 4 files changed, 46 insertions(+), 22 deletions(-) diff --git a/src/lp_data/HStruct.h b/src/lp_data/HStruct.h index 38e78855df..d69c36dede 100644 --- a/src/lp_data/HStruct.h +++ b/src/lp_data/HStruct.h @@ -58,6 +58,16 @@ struct HotStart { }; struct HighsBasis { + // Logical flags for a HiGHS basis: + // + // valid: has been factored by HiGHS + // + // alien: a basis that's been set externally, so cannot be assumed + // to even have the right number of basic and nonbasic variables + // + // useful: a basis that may be useful + // + // Need useful since, by default, a basis is alien but not useful bool valid = false; bool alien = true; bool useful = false; @@ -67,8 +77,8 @@ struct HighsBasis { std::string debug_origin_name = "None"; std::vector col_status; std::vector row_status; - void print() const; - void printScalars() const; + void print(std::string message = "") const; + void printScalars(std::string message = "") const; void invalidate(); void clear(); }; diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index 3e45df3443..6ad3dd9c77 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -1266,19 +1266,14 @@ HighsStatus Highs::solve() { // is available, simplex should surely be chosen. const bool solver_will_use_basis = options_.solver == kSimplexString || options_.solver == kHighsChooseString; - const bool has_basis = basis_.valid && basis_.useful; + const bool has_basis = basis_.valid || basis_.useful; if (has_basis) { assert(basis_.col_status.size() == static_cast(incumbent_lp.num_col_)); assert(basis_.row_status.size() == static_cast(incumbent_lp.num_row_)); } if (basis_.valid) assert(basis_.useful); - if (has_basis && !basis_.valid) { - basis_.print(); - printf("Highs::solve() has_basis && !basis_.valid\n"); - assert(111==123); - } - if ((basis_.valid || options_.presolve == kHighsOffString || + if ((has_basis || options_.presolve == kHighsOffString || unconstrained_lp) && solver_will_use_basis) { // There is a valid basis for the problem, presolve is off, or LP diff --git a/src/lp_data/HighsSolution.cpp b/src/lp_data/HighsSolution.cpp index fb08dc9c2f..fe45531e52 100644 --- a/src/lp_data/HighsSolution.cpp +++ b/src/lp_data/HighsSolution.cpp @@ -1357,9 +1357,15 @@ HighsStatus formSimplexLpBasisAndFactor(HighsLpSolverObject& solver_object, // If new scaling is performed, the hot start information is // no longer valid if (new_scaling) ekk_instance.clearHotStart(); - if (basis.alien) { - // An alien basis needs to be checked for rank deficiency, and - // possibly completed if it is rectangular + const bool check_basis = basis.alien || + (!basis.valid && basis.useful); + if (check_basis) { + // The basis needs to be checked for rank deficiency, and possibly + // completed if it is rectangular + // + // If it's not valid but useful, but not alien, + // accommodateAlienBasis will assert, so make the basis alien + basis.alien = true; assert(!only_from_known_basis); accommodateAlienBasis(solver_object); basis.alien = false; @@ -1606,23 +1612,24 @@ void HighsSolution::clear() { void HighsObjectiveSolution::clear() { this->col_value.clear(); } -void HighsBasis::print() const { +void HighsBasis::print(std::string message) const { if (!this->useful) return; - this->printScalars(); + this->printScalars(message); for (HighsInt iCol = 0; iCol < HighsInt(this->col_status.size()); iCol++) printf("Basis: col_status[%2d] = %d\n", int(iCol), int(this->col_status[iCol])); for (HighsInt iRow = 0; iRow < HighsInt(this->row_status.size()); iRow++) printf("Basis: row_status[%2d] = %d\n", int(iRow), int(this->row_status[iRow])); } -void HighsBasis::printScalars() const { - printf("Basis: valid = %d\n", this->valid); - printf("Basis: alien = %d\n", this->alien); - printf("Basis: useful = %d\n", this->useful); - printf("Basis: was_alien = %d\n", this->was_alien); - printf("Basis: debug_id = %d\n", int(this->debug_id)); - printf("Basis: debug_update_count = %d\n", int(this->debug_update_count)); - printf("Basis: debug_origin_name = %s\n", this->debug_origin_name.c_str()); +void HighsBasis::printScalars(std::string message) const { + printf("\nBasis: %s\n", message.c_str()); + printf(" valid = %d\n", this->valid); + printf(" alien = %d\n", this->alien); + printf(" useful = %d\n", this->useful); + printf(" was_alien = %d\n", this->was_alien); + printf(" debug_id = %d\n", int(this->debug_id)); + printf(" debug_update_count = %d\n", int(this->debug_update_count)); + printf(" debug_origin_name = %s\n", this->debug_origin_name.c_str()); } void HighsBasis::invalidate() { diff --git a/src/simplex/HApp.h b/src/simplex/HApp.h index 9fc35589b4..a19d6080e0 100644 --- a/src/simplex/HApp.h +++ b/src/simplex/HApp.h @@ -136,6 +136,18 @@ inline HighsStatus solveLpSimplex(HighsLpSolverObject& solver_object) { // If new scaling is performed, the hot start information is // no longer valid if (new_scaling) ekk_instance.clearHotStart(); + // + if (!status.has_basis && !basis.valid && basis.useful) { + // There is no simplex basis, but there is a useful HiGHS basis + // that is not validated + assert(basis.col_status.size() == static_cast(incumbent_lp.num_col_)); + assert(basis.row_status.size() == static_cast(incumbent_lp.num_row_)); + basis.print("Before formSimplexLpBasisAndFactor"); + HighsStatus return_status = formSimplexLpBasisAndFactor(solver_object); + if (return_status != HighsStatus::kOk) return returnFromSolveLpSimplex(solver_object, HighsStatus::kError); + basis.print("After formSimplexLpBasisAndFactor"); + basis.valid = true; + } // Move the LP to EKK, updating other EKK pointers and any simplex // NLA pointers, since they may have moved if the LP has been // modified From 4f21f0a734cfbb7348d46b808ad5889b701ae778 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 18:20:22 +0000 Subject: [PATCH 06/11] No longer necessarily invalidating basis in deleteCols; all unit tests pass --- src/lp_data/Highs.cpp | 2 +- src/lp_data/HighsInterface.cpp | 14 +++++++++----- 2 files changed, 10 insertions(+), 6 deletions(-) diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index 6ad3dd9c77..4597573bbb 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -1266,7 +1266,7 @@ HighsStatus Highs::solve() { // is available, simplex should surely be chosen. const bool solver_will_use_basis = options_.solver == kSimplexString || options_.solver == kHighsChooseString; - const bool has_basis = basis_.valid || basis_.useful; + const bool has_basis = basis_.useful; if (has_basis) { assert(basis_.col_status.size() == static_cast(incumbent_lp.num_col_)); assert(basis_.row_status.size() == static_cast(incumbent_lp.num_row_)); diff --git a/src/lp_data/HighsInterface.cpp b/src/lp_data/HighsInterface.cpp index 31eb3edcd2..54fa94eacc 100644 --- a/src/lp_data/HighsInterface.cpp +++ b/src/lp_data/HighsInterface.cpp @@ -601,8 +601,6 @@ void deleteBasisEntries(std::vector& status, keep_from_entry, keep_to_entry, current_set_entry); // Account for the initial entries being kept if (k == from_k) new_num_entry = delete_from_entry; - if (delete_to_entry >= entry_dim - 1) break; - assert(delete_to_entry < entry_dim); // Identify whether a basic or a nonbasic entry has been deleted for (HighsInt entry = delete_from_entry; entry <= delete_to_entry; entry++) { if (status[entry] == HighsBasisStatus::kBasic) { @@ -611,6 +609,7 @@ void deleteBasisEntries(std::vector& status, deleted_nonbasic = true; } } + if (delete_to_entry >= entry_dim - 1) break; for (HighsInt entry = keep_from_entry; entry <= keep_to_entry; entry++) { status[new_num_entry] = status[entry]; new_num_entry++; @@ -674,8 +673,7 @@ void Highs::deleteColsInterface(HighsIndexCollection& index_collection) { } else { assert(!basis.valid); } - basis.valid = false; - + if (lp.scale_.has_scaling) { deleteScale(lp.scale_.col, index_collection); lp.scale_.col.resize(lp.num_col_); @@ -2057,7 +2055,13 @@ HighsStatus Highs::elasticityFilterReturn( run_status = this->deleteCols(original_num_col, lp.num_col_ - 1); assert(run_status == HighsStatus::kOk); - + // + // Now that deleteRows and deleteCols may yield a valid basis, the + // lack of dual values triggers an assert in + // getKktFailures. Ultimately (#2081) the dual values will be + // available but, for now, make the basis invalid. + basis_.valid = false; + run_status = this->changeColsCost(0, original_num_col - 1, original_col_cost.data()); assert(run_status == HighsStatus::kOk); From 20fb52db768df931dc3739e9cc93271604d0b23c Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 18:31:16 +0000 Subject: [PATCH 07/11] No longer necessarily invalidating basis in deleteRows; all unit tests pass --- src/lp_data/HighsInterface.cpp | 1 - src/simplex/HApp.h | 2 -- 2 files changed, 3 deletions(-) diff --git a/src/lp_data/HighsInterface.cpp b/src/lp_data/HighsInterface.cpp index 54fa94eacc..635298e882 100644 --- a/src/lp_data/HighsInterface.cpp +++ b/src/lp_data/HighsInterface.cpp @@ -729,7 +729,6 @@ void Highs::deleteRowsInterface(HighsIndexCollection& index_collection) { } else { assert(!basis.valid); } - basis.valid = false; if (lp.scale_.has_scaling) { deleteScale(lp.scale_.row, index_collection); diff --git a/src/simplex/HApp.h b/src/simplex/HApp.h index a19d6080e0..b7d073bdea 100644 --- a/src/simplex/HApp.h +++ b/src/simplex/HApp.h @@ -142,10 +142,8 @@ inline HighsStatus solveLpSimplex(HighsLpSolverObject& solver_object) { // that is not validated assert(basis.col_status.size() == static_cast(incumbent_lp.num_col_)); assert(basis.row_status.size() == static_cast(incumbent_lp.num_row_)); - basis.print("Before formSimplexLpBasisAndFactor"); HighsStatus return_status = formSimplexLpBasisAndFactor(solver_object); if (return_status != HighsStatus::kOk) return returnFromSolveLpSimplex(solver_object, HighsStatus::kError); - basis.print("After formSimplexLpBasisAndFactor"); basis.valid = true; } // Move the LP to EKK, updating other EKK pointers and any simplex From 78793f96a4588c57948e49fddb82c4169f604302 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 18:55:10 +0000 Subject: [PATCH 08/11] Now taking 1 iteration in bin/unit_tests hot-start-after-delete for nonbasic row deletion and addition --- src/lp_data/Highs.cpp | 2 +- src/lp_data/HighsInterface.cpp | 24 +++++++++--------------- src/lp_data/HighsSolution.cpp | 2 +- src/simplex/HApp.h | 3 +++ 4 files changed, 14 insertions(+), 17 deletions(-) diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index 4597573bbb..eace21d527 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -1282,7 +1282,7 @@ HighsStatus Highs::solve() { "LP without presolve, or with basis, or unconstrained"; // If there is a valid HiGHS basis, refine any status values that // are simply HighsBasisStatus::kNonbasic - if (basis_.valid) refineBasis(incumbent_lp, solution_, basis_); + if (basis_.useful) refineBasis(incumbent_lp, solution_, basis_); solveLp(incumbent_lp, "Solving LP without presolve, or with basis, or unconstrained", this_solve_original_lp_time); diff --git a/src/lp_data/HighsInterface.cpp b/src/lp_data/HighsInterface.cpp index 635298e882..e5a1c50fc7 100644 --- a/src/lp_data/HighsInterface.cpp +++ b/src/lp_data/HighsInterface.cpp @@ -471,7 +471,8 @@ HighsStatus Highs::addRowsInterface(HighsInt ext_num_new_row, HighsLp& lp = model_.lp_; HighsBasis& basis = basis_; HighsScale& scale = lp.scale_; - bool& valid_basis = basis.valid; + bool& useful_basis = basis.useful; + bool& lp_has_scaling = lp.scale_.has_scaling; // Check that if nonzeros are to be added then the model has a positive number @@ -559,7 +560,7 @@ HighsStatus Highs::addRowsInterface(HighsInt ext_num_new_row, &scale.row[lp.num_row_]); } // Update the basis corresponding to new basic rows - if (valid_basis) appendBasicRowsToBasisInterface(ext_num_new_row); + if (useful_basis) appendBasicRowsToBasisInterface(ext_num_new_row); // Possibly add row names lp.addRowNames("", ext_num_new_row); @@ -1343,27 +1344,20 @@ void Highs::appendNonbasicColsToBasisInterface(const HighsInt ext_num_new_col) { void Highs::appendBasicRowsToBasisInterface(const HighsInt ext_num_new_row) { if (ext_num_new_row == 0) return; HighsBasis& highs_basis = basis_; - if (!highs_basis.valid) return; + if (!highs_basis.useful) return; const bool has_simplex_basis = ekk_instance_.status_.has_basis; SimplexBasis& simplex_basis = ekk_instance_.basis_; HighsLp& lp = model_.lp_; - const bool has_highs_basis = highs_basis.valid && highs_basis.useful; - - if (has_highs_basis) { - assert(highs_basis.col_status.size() == static_cast(lp.num_col_)); - assert(highs_basis.row_status.size() == static_cast(lp.num_row_)); - } - if (!has_highs_basis && !has_simplex_basis) return; + assert(highs_basis.col_status.size() == static_cast(lp.num_col_)); + assert(highs_basis.row_status.size() == static_cast(lp.num_row_)); // Add basic logicals // Add the new rows to the Highs basis HighsInt newNumRow = lp.num_row_ + ext_num_new_row; - if (has_highs_basis) { - highs_basis.row_status.resize(newNumRow); - for (HighsInt iRow = lp.num_row_; iRow < newNumRow; iRow++) - highs_basis.row_status[iRow] = HighsBasisStatus::kBasic; - } + highs_basis.row_status.resize(newNumRow); + for (HighsInt iRow = lp.num_row_; iRow < newNumRow; iRow++) + highs_basis.row_status[iRow] = HighsBasisStatus::kBasic; if (has_simplex_basis) { // Add the new rows to the simplex basis HighsInt newNumTot = lp.num_col_ + newNumRow; diff --git a/src/lp_data/HighsSolution.cpp b/src/lp_data/HighsSolution.cpp index fe45531e52..7dee62abe5 100644 --- a/src/lp_data/HighsSolution.cpp +++ b/src/lp_data/HighsSolution.cpp @@ -677,7 +677,7 @@ double computeObjectiveValue(const HighsLp& lp, const HighsSolution& solution) { // and any solution values void refineBasis(const HighsLp& lp, const HighsSolution& solution, HighsBasis& basis) { - assert(basis.valid); + assert(basis.useful); assert(isBasisRightSize(lp, basis)); const bool have_highs_solution = solution.value_valid; diff --git a/src/simplex/HApp.h b/src/simplex/HApp.h index b7d073bdea..8b63ce74b8 100644 --- a/src/simplex/HApp.h +++ b/src/simplex/HApp.h @@ -144,6 +144,9 @@ inline HighsStatus solveLpSimplex(HighsLpSolverObject& solver_object) { assert(basis.row_status.size() == static_cast(incumbent_lp.num_row_)); HighsStatus return_status = formSimplexLpBasisAndFactor(solver_object); if (return_status != HighsStatus::kOk) return returnFromSolveLpSimplex(solver_object, HighsStatus::kError); + // formSimplexLpBasisAndFactor may introduce variables with + // HighsBasisStatus::kNonbasic, so refine it + refineBasis(incumbent_lp, solution, basis); basis.valid = true; } // Move the LP to EKK, updating other EKK pointers and any simplex From 4238871f1aa8bd9b3717131281a9225ca6c32055 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 19:04:16 +0000 Subject: [PATCH 09/11] Now calling appendNonbasicColsToBasisInterface when basis is useful, and this only assumes a useful basis --- src/lp_data/HighsInterface.cpp | 18 +++++++----------- 1 file changed, 7 insertions(+), 11 deletions(-) diff --git a/src/lp_data/HighsInterface.cpp b/src/lp_data/HighsInterface.cpp index e5a1c50fc7..e7faea3860 100644 --- a/src/lp_data/HighsInterface.cpp +++ b/src/lp_data/HighsInterface.cpp @@ -305,7 +305,7 @@ HighsStatus Highs::addColsInterface( HighsLp& lp = model_.lp_; HighsBasis& basis = basis_; HighsScale& scale = lp.scale_; - bool& valid_basis = basis.valid; + bool& useful_basis = basis.useful; bool& lp_has_scaling = lp.scale_.has_scaling; // Check that if nonzeros are to be added then the model has a positive number @@ -415,7 +415,7 @@ HighsStatus Highs::addColsInterface( &scale.col[lp.num_col_]); } // Update the basis corresponding to new nonbasic columns - if (valid_basis) appendNonbasicColsToBasisInterface(ext_num_new_col); + if (useful_basis) appendNonbasicColsToBasisInterface(ext_num_new_col); // Possibly add column names lp.addColNames("", ext_num_new_col); @@ -1263,22 +1263,18 @@ void Highs::setNonbasicStatusInterface( void Highs::appendNonbasicColsToBasisInterface(const HighsInt ext_num_new_col) { if (ext_num_new_col == 0) return; HighsBasis& highs_basis = basis_; - if (!highs_basis.valid) return; + if (!highs_basis.useful) return; const bool has_simplex_basis = ekk_instance_.status_.has_basis; SimplexBasis& simplex_basis = ekk_instance_.basis_; HighsLp& lp = model_.lp_; - const bool has_highs_basis = highs_basis.valid && highs_basis.useful; - if (has_highs_basis) { - assert(highs_basis.col_status.size() == static_cast(lp.num_col_)); - assert(highs_basis.row_status.size() == static_cast(lp.num_row_)); - } - if (!has_highs_basis && !has_simplex_basis) return; + assert(highs_basis.col_status.size() == static_cast(lp.num_col_)); + assert(highs_basis.row_status.size() == static_cast(lp.num_row_)); // Add nonbasic structurals HighsInt newNumCol = lp.num_col_ + ext_num_new_col; HighsInt newNumTot = newNumCol + lp.num_row_; - if (has_highs_basis) highs_basis.col_status.resize(newNumCol); + highs_basis.col_status.resize(newNumCol); if (has_simplex_basis) { simplex_basis.nonbasicFlag_.resize(newNumTot); simplex_basis.nonbasicMove_.resize(newNumTot); @@ -1332,7 +1328,7 @@ void Highs::appendNonbasicColsToBasisInterface(const HighsInt ext_num_new_col) { move = kNonbasicMoveZe; } assert(status != HighsBasisStatus::kNonbasic); - if (has_highs_basis) highs_basis.col_status[iCol] = status; + highs_basis.col_status[iCol] = status; if (has_simplex_basis) { assert(move != kIllegalMoveValue); simplex_basis.nonbasicFlag_[iCol] = kNonbasicFlagTrue; From 407979409f6a147691d3f553bc7fc3e6fab0de1a Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 21:44:39 +0000 Subject: [PATCH 10/11] Added testing with set and mask; formatted --- check/TestLpModification.cpp | 223 ++++++++++++++++++++++++++------- src/lp_data/Highs.cpp | 9 +- src/lp_data/HighsInterface.cpp | 52 ++++---- src/lp_data/HighsSolution.cpp | 9 +- src/simplex/HApp.h | 11 +- 5 files changed, 221 insertions(+), 83 deletions(-) diff --git a/check/TestLpModification.cpp b/check/TestLpModification.cpp index be24bf6605..9a443ea43f 100644 --- a/check/TestLpModification.cpp +++ b/check/TestLpModification.cpp @@ -2273,7 +2273,6 @@ TEST_CASE("row-wise-get-row-avgas", "[highs_data]") { h.ensureColwise(); testAvgasGetRow(h); testAvgasGetCol(h); - } TEST_CASE("hot-start-after-delete", "[highs_data]") { @@ -2285,54 +2284,194 @@ TEST_CASE("hot-start-after-delete", "[highs_data]") { const HighsSolution& solution = h.getSolution(); std::string model = "avgas"; std::string model_file = - std::string(HIGHS_DIR) + "/check/instances/" + model + ".mps"; + std::string(HIGHS_DIR) + "/check/instances/" + model + ".mps"; h.readModel(model_file); h.run(); - printf("Initial solve takes %d iterations and yields objective = %g\n", - int(info.simplex_iteration_count), info.objective_function_value); - h.writeSolution("", 1); - for (HighsInt iCol = 0; iCol < lp.num_col_; iCol++) - printf("Col %2d is %s\n", int(iCol), basis.col_status[iCol] == HighsBasisStatus::kBasic ? "basic" : "nonbasic"); - printf("\n"); - for (HighsInt iRow = 0; iRow < lp.num_row_; iRow++) - printf("Row %2d is %s\n", int(iRow), basis.row_status[iRow] == HighsBasisStatus::kBasic ? "basic" : "nonbasic"); + HighsInt ieration_count0 = info.simplex_iteration_count; + if (dev_run) + printf("Initial solve takes %d iterations and yields objective = %g\n", + int(info.simplex_iteration_count), info.objective_function_value); - double cost; - double lower; - double upper; + HighsInt max_dim = std::max(lp.num_col_, lp.num_row_); + std::vector cost(1); + std::vector lower(1); + std::vector upper(1); HighsInt nnz; std::vector start(1); - std::vector index(lp.num_row_); - std::vector value(lp.num_row_); + std::vector index(max_dim); + std::vector value(max_dim); HighsInt get_num; + HighsInt use_col, use_row; + for (HighsInt k = 0; k < 2; k++) { + if (dev_run) { + for (HighsInt iCol = 0; iCol < lp.num_col_; iCol++) + printf("Col %2d is %s\n", int(iCol), + basis.col_status[iCol] == HighsBasisStatus::kBasic ? "basic" + : "nonbasic"); + printf("\n"); + } + if (k == 0) { + use_col = 1; // Nonbasic + } else { + use_col = 4; // Basic + } + if (dev_run) + printf( + "\nDeleting and adding column %1d with status \"%s\" and value %g\n", + int(use_col), + h.basisStatusToString(basis.col_status[use_col]).c_str(), + solution.col_value[use_col]); + + h.getCols(use_col, use_col, get_num, cost.data(), lower.data(), + upper.data(), nnz, start.data(), index.data(), value.data()); + + h.deleteCols(use_col, use_col); + if (dev_run) basis.printScalars(); + + h.addCol(cost[0], lower[0], upper[0], nnz, index.data(), value.data()); + + h.run(); + if (dev_run) + printf( + "After deleting and adding column %1d, solve takes %d iterations and " + "yields objective = %g\n", + int(use_col), int(info.simplex_iteration_count), + info.objective_function_value); + REQUIRE(info.simplex_iteration_count < ieration_count0); + } + + for (HighsInt k = 0; k < 2; k++) { + if (dev_run) { + for (HighsInt iRow = 0; iRow < lp.num_row_; iRow++) + printf("Row %2d is %s\n", int(iRow), + basis.row_status[iRow] == HighsBasisStatus::kBasic ? "basic" + : "nonbasic"); + } + if (k == 0) { + use_row = 1; // Nonbasic + } else { + use_row = 8; // Basic + } + if (dev_run) + printf("\nDeleting and adding row %1d with status \"%s\" and value %g\n", + int(use_row), + h.basisStatusToString(basis.row_status[use_row]).c_str(), + solution.row_value[use_row]); + + h.getRows(use_row, use_row, get_num, lower.data(), upper.data(), nnz, + start.data(), index.data(), value.data()); + + h.deleteRows(use_row, use_row); + if (dev_run) basis.printScalars(); + + h.addRow(lower[0], upper[0], nnz, index.data(), value.data()); + + h.run(); + if (dev_run) + printf( + "After deleting and adding row %1d, solve takes %d iterations and " + "yields objective = %g\n", + int(use_row), int(info.simplex_iteration_count), + info.objective_function_value); + REQUIRE(info.simplex_iteration_count < ieration_count0); + } + std::vector set = {1, 3, 4}; + HighsInt num_set_en = set.size(); + cost.resize(num_set_en); + lower.resize(num_set_en); + upper.resize(num_set_en); + start.resize(num_set_en); + index.resize(num_set_en * max_dim); + value.resize(num_set_en * max_dim); + + h.getCols(num_set_en, set.data(), get_num, cost.data(), lower.data(), + upper.data(), nnz, start.data(), index.data(), value.data()); + + h.deleteCols(num_set_en, set.data()); + if (dev_run) basis.printScalars(); + + h.addCols(get_num, cost.data(), lower.data(), upper.data(), nnz, start.data(), + index.data(), value.data()); + + h.run(); + if (dev_run) + printf( + "After deleting and adding %d columns in set, solve takes %d " + "iterations and yields objective = %g\n", + int(get_num), int(info.simplex_iteration_count), + info.objective_function_value); + // REQUIRE(info.simplex_iteration_count < ieration_count0); + + h.getRows(num_set_en, set.data(), get_num, lower.data(), upper.data(), nnz, + start.data(), index.data(), value.data()); + + h.deleteRows(num_set_en, set.data()); + if (dev_run) basis.printScalars(); + + h.addRows(get_num, lower.data(), upper.data(), nnz, start.data(), + index.data(), value.data()); + + h.run(); + if (dev_run) + printf( + "After deleting and adding %d rows in set, solve takes %d iterations " + "and yields objective = %g\n", + int(get_num), int(info.simplex_iteration_count), + info.objective_function_value); + // REQUIRE(info.simplex_iteration_count < ieration_count0); + std::vector mask; + mask.assign(max_dim, 0); + mask[1] = 1; + mask[4] = 1; + mask[5] = 1; + + h.getCols(mask.data(), get_num, cost.data(), lower.data(), upper.data(), nnz, + start.data(), index.data(), value.data()); + + h.deleteCols(mask.data()); + if (dev_run) basis.printScalars(); + + h.addCols(get_num, cost.data(), lower.data(), upper.data(), nnz, start.data(), + index.data(), value.data()); - HighsInt use_col = 1; - printf("\nDeleting and adding column %1d with status \"%s\" and value %g\n", - int(use_col), h.basisStatusToString(basis.col_status[use_col]).c_str(), solution.col_value[use_col]); - - h.getCols(use_col, use_col, get_num, &cost, &lower, &upper, nnz, start.data(), index.data(), value.data()); - - h.deleteCols(use_col, use_col); - basis.printScalars(); - - h.addCol(cost, lower, upper, nnz, index.data(), value.data()); - h.run(); - printf("After deleting and adding column %1d, solve takes %d iterations and yields objective = %g\n", - int(use_col), int(info.simplex_iteration_count), info.objective_function_value); - - HighsInt use_row = 1; - printf("\nDeleting and adding row %1d with status \"%s\" and value %g\n", - int(use_row), h.basisStatusToString(basis.row_status[use_row]).c_str(), solution.row_value[use_row]); - - h.getRows(use_row, use_row, get_num, &lower, &upper, nnz, start.data(), index.data(), value.data()); - - h.deleteRows(use_row, use_row); - basis.printScalars(); - - h.addRow(lower, upper, nnz, index.data(), value.data()); - + if (dev_run) + printf( + "After deleting and adding %d columns in mask, solve takes %d " + "iterations and yields objective = %g\n", + int(get_num), int(info.simplex_iteration_count), + info.objective_function_value); + // REQUIRE(info.simplex_iteration_count < ieration_count0); + + mask.assign(max_dim, 0); + mask[1] = 1; + mask[4] = 1; + mask[5] = 1; + mask[8] = 1; + mask[9] = 1; + HighsInt num_mask_en = mask.size(); + cost.resize(num_mask_en); + lower.resize(num_mask_en); + upper.resize(num_mask_en); + start.resize(num_mask_en); + index.resize(num_mask_en * max_dim); + value.resize(num_mask_en * max_dim); + + h.getRows(mask.data(), get_num, lower.data(), upper.data(), nnz, start.data(), + index.data(), value.data()); + + h.deleteRows(mask.data()); + if (dev_run) basis.printScalars(); + + h.addRows(get_num, lower.data(), upper.data(), nnz, start.data(), + index.data(), value.data()); + h.run(); - printf("After deleting and adding row %1d, solve takes %d iterations and yields objective = %g\n", - int(use_row), int(info.simplex_iteration_count), info.objective_function_value); + if (dev_run) + printf( + "After deleting and adding %d rows in mask, solve takes %d iterations " + "and yields objective = %g\n", + int(get_num), int(info.simplex_iteration_count), + info.objective_function_value); + // REQUIRE(info.simplex_iteration_count < ieration_count0); } diff --git a/src/lp_data/Highs.cpp b/src/lp_data/Highs.cpp index eace21d527..8372d4aded 100644 --- a/src/lp_data/Highs.cpp +++ b/src/lp_data/Highs.cpp @@ -1268,13 +1268,14 @@ HighsStatus Highs::solve() { options_.solver == kHighsChooseString; const bool has_basis = basis_.useful; if (has_basis) { - assert(basis_.col_status.size() == static_cast(incumbent_lp.num_col_)); - assert(basis_.row_status.size() == static_cast(incumbent_lp.num_row_)); + assert(basis_.col_status.size() == + static_cast(incumbent_lp.num_col_)); + assert(basis_.row_status.size() == + static_cast(incumbent_lp.num_row_)); } if (basis_.valid) assert(basis_.useful); - if ((has_basis || options_.presolve == kHighsOffString || - unconstrained_lp) && + if ((has_basis || options_.presolve == kHighsOffString || unconstrained_lp) && solver_will_use_basis) { // There is a valid basis for the problem, presolve is off, or LP // has no constraint matrix, and the solver will use the basis diff --git a/src/lp_data/HighsInterface.cpp b/src/lp_data/HighsInterface.cpp index e7faea3860..627338c7e1 100644 --- a/src/lp_data/HighsInterface.cpp +++ b/src/lp_data/HighsInterface.cpp @@ -472,7 +472,7 @@ HighsStatus Highs::addRowsInterface(HighsInt ext_num_new_row, HighsBasis& basis = basis_; HighsScale& scale = lp.scale_; bool& useful_basis = basis.useful; - + bool& lp_has_scaling = lp.scale_.has_scaling; // Check that if nonzeros are to be added then the model has a positive number @@ -578,10 +578,9 @@ HighsStatus Highs::addRowsInterface(HighsInt ext_num_new_row, } void deleteBasisEntries(std::vector& status, - bool& deleted_basic, - bool& deleted_nonbasic, - const HighsIndexCollection& index_collection, - const HighsInt entry_dim) { + bool& deleted_basic, bool& deleted_nonbasic, + const HighsIndexCollection& index_collection, + const HighsInt entry_dim) { assert(ok(index_collection)); assert(static_cast(entry_dim) == status.size()); HighsInt from_k; @@ -603,11 +602,12 @@ void deleteBasisEntries(std::vector& status, // Account for the initial entries being kept if (k == from_k) new_num_entry = delete_from_entry; // Identify whether a basic or a nonbasic entry has been deleted - for (HighsInt entry = delete_from_entry; entry <= delete_to_entry; entry++) { + for (HighsInt entry = delete_from_entry; entry <= delete_to_entry; + entry++) { if (status[entry] == HighsBasisStatus::kBasic) { - deleted_basic = true; + deleted_basic = true; } else { - deleted_nonbasic = true; + deleted_nonbasic = true; } } if (delete_to_entry >= entry_dim - 1) break; @@ -620,30 +620,24 @@ void deleteBasisEntries(std::vector& status, status.resize(new_num_entry); } -void deleteBasisCols(HighsBasis& basis, - const HighsIndexCollection& index_collection, - const HighsInt original_num_col) { +void deleteBasisCols(HighsBasis& basis, + const HighsIndexCollection& index_collection, + const HighsInt original_num_col) { bool deleted_basic; bool deleted_nonbasic; - deleteBasisEntries(basis.col_status, - deleted_basic, - deleted_nonbasic, - index_collection, original_num_col); - if (deleted_basic) - basis.valid = false; + deleteBasisEntries(basis.col_status, deleted_basic, deleted_nonbasic, + index_collection, original_num_col); + if (deleted_basic) basis.valid = false; } -void deleteBasisRows(HighsBasis& basis, - const HighsIndexCollection& index_collection, - const HighsInt original_num_row) { +void deleteBasisRows(HighsBasis& basis, + const HighsIndexCollection& index_collection, + const HighsInt original_num_row) { bool deleted_basic; bool deleted_nonbasic; - deleteBasisEntries(basis.row_status, - deleted_basic, - deleted_nonbasic, - index_collection, original_num_row); - if (deleted_nonbasic) - basis.valid = false; + deleteBasisEntries(basis.row_status, deleted_basic, deleted_nonbasic, + index_collection, original_num_row); + if (deleted_nonbasic) basis.valid = false; } void Highs::deleteColsInterface(HighsIndexCollection& index_collection) { @@ -674,7 +668,7 @@ void Highs::deleteColsInterface(HighsIndexCollection& index_collection) { } else { assert(!basis.valid); } - + if (lp.scale_.has_scaling) { deleteScale(lp.scale_.col, index_collection); lp.scale_.col.resize(lp.num_col_); @@ -2044,13 +2038,13 @@ HighsStatus Highs::elasticityFilterReturn( run_status = this->deleteCols(original_num_col, lp.num_col_ - 1); assert(run_status == HighsStatus::kOk); - // + // // Now that deleteRows and deleteCols may yield a valid basis, the // lack of dual values triggers an assert in // getKktFailures. Ultimately (#2081) the dual values will be // available but, for now, make the basis invalid. basis_.valid = false; - + run_status = this->changeColsCost(0, original_num_col - 1, original_col_cost.data()); assert(run_status == HighsStatus::kOk); diff --git a/src/lp_data/HighsSolution.cpp b/src/lp_data/HighsSolution.cpp index 7dee62abe5..a28060861a 100644 --- a/src/lp_data/HighsSolution.cpp +++ b/src/lp_data/HighsSolution.cpp @@ -1357,8 +1357,7 @@ HighsStatus formSimplexLpBasisAndFactor(HighsLpSolverObject& solver_object, // If new scaling is performed, the hot start information is // no longer valid if (new_scaling) ekk_instance.clearHotStart(); - const bool check_basis = basis.alien || - (!basis.valid && basis.useful); + const bool check_basis = basis.alien || (!basis.valid && basis.useful); if (check_basis) { // The basis needs to be checked for rank deficiency, and possibly // completed if it is rectangular @@ -1616,9 +1615,11 @@ void HighsBasis::print(std::string message) const { if (!this->useful) return; this->printScalars(message); for (HighsInt iCol = 0; iCol < HighsInt(this->col_status.size()); iCol++) - printf("Basis: col_status[%2d] = %d\n", int(iCol), int(this->col_status[iCol])); + printf("Basis: col_status[%2d] = %d\n", int(iCol), + int(this->col_status[iCol])); for (HighsInt iRow = 0; iRow < HighsInt(this->row_status.size()); iRow++) - printf("Basis: row_status[%2d] = %d\n", int(iRow), int(this->row_status[iRow])); + printf("Basis: row_status[%2d] = %d\n", int(iRow), + int(this->row_status[iRow])); } void HighsBasis::printScalars(std::string message) const { diff --git a/src/simplex/HApp.h b/src/simplex/HApp.h index 8b63ce74b8..178593cbc6 100644 --- a/src/simplex/HApp.h +++ b/src/simplex/HApp.h @@ -136,14 +136,17 @@ inline HighsStatus solveLpSimplex(HighsLpSolverObject& solver_object) { // If new scaling is performed, the hot start information is // no longer valid if (new_scaling) ekk_instance.clearHotStart(); - // + // if (!status.has_basis && !basis.valid && basis.useful) { // There is no simplex basis, but there is a useful HiGHS basis // that is not validated - assert(basis.col_status.size() == static_cast(incumbent_lp.num_col_)); - assert(basis.row_status.size() == static_cast(incumbent_lp.num_row_)); + assert(basis.col_status.size() == + static_cast(incumbent_lp.num_col_)); + assert(basis.row_status.size() == + static_cast(incumbent_lp.num_row_)); HighsStatus return_status = formSimplexLpBasisAndFactor(solver_object); - if (return_status != HighsStatus::kOk) return returnFromSolveLpSimplex(solver_object, HighsStatus::kError); + if (return_status != HighsStatus::kOk) + return returnFromSolveLpSimplex(solver_object, HighsStatus::kError); // formSimplexLpBasisAndFactor may introduce variables with // HighsBasisStatus::kNonbasic, so refine it refineBasis(incumbent_lp, solution, basis); From f3664d9e023ba3981fe3c314ddac374e7adf0056 Mon Sep 17 00:00:00 2001 From: JAJHall Date: Mon, 13 Jan 2025 22:17:55 +0000 Subject: [PATCH 11/11] Updated FEATURES.md --- FEATURES.md | 2 ++ 1 file changed, 2 insertions(+) diff --git a/FEATURES.md b/FEATURES.md index eb17f2ed97..f6ff74c98c 100644 --- a/FEATURES.md +++ b/FEATURES.md @@ -17,3 +17,5 @@ Fixed bug in presolve when pointers stored in HighsMatrixSlice get invalidated w Primal and dual residual tolerances - applied following IPM or PDLP solution - now documented as options Highs::getCols (Highs::getRows) now runs in linear time if the internal constraint matrix is stored column-wise (row-wise). Added ensureColwise/Rowwise to the Highs class, the C API and highspy so that users can set the internal constraint matrix storage orientation + +When columns and rows are deleted from the incumbent LP after a basic solution has been found, HiGHS no longer invalidates the basis. Now it maintains the basic and nonbasic status of the remaining variables and constraints. When the model is re-solved, this information is used to construct a starting basis. \ No newline at end of file