Need HighsIis::feasible_

Author: jajhall

check/TestIis.cpp

@@ -491,7 +491,6 @@ TEST_CASE("lp-feasibility-relaxation", "[iis]") {
   h.setOptionValue("output_flag", dev_run);
   const HighsSolution& solution = h.getSolution();
   h.passModel(lp);
-  //  h.run();
 
   const bool all_tests = false;
   const bool test0 = false || all_tests;
@@ -586,7 +585,7 @@ void testMps(std::string& model, const HighsInt iis_strategy,
   std::string model_file =
       std::string(HIGHS_DIR) + "/check/instances/" + model + ".mps";
   Highs highs;
-  //  highs.setOptionValue("output_flag", dev_run);
+  highs.setOptionValue("output_flag", dev_run);
 
   REQUIRE(highs.readModel(model_file) == HighsStatus::kOk);
   //  if (iis_strategy == kIisStrategyFromRay ||
@@ -658,7 +657,7 @@ TEST_CASE("feasible-lp-iis", "[iis]") {
   lp.a_matrix_.index_ = {0, 1, 0, 1};
   lp.a_matrix_.value_ = {1, 1, 2, 4};
   Highs h;
-  h.setOptionValue("output_flag", dev_run);
+  //  h.setOptionValue("output_flag", dev_run);
   h.passModel(lp);
   HighsIis iis;
   // With kIisStrategyLight, feasibility of the LP is not determined

highs/lp_data/HighsIis.cpp

@@ -469,10 +469,28 @@ void HighsIis::setLp(const HighsLp& lp) {
   iis_lp.model_name_ = lp.model_name_ + "_IIS";
 }
 
-void HighsIis::setStatus(const HighsLp& lp) {
+HighsInt HighsIis::nonIsStatus(const HighsModelStatus& model_status) const {
+  const bool is_feasible =
+    model_status == HighsModelStatus::kOptimal ||
+    model_status == HighsModelStatus::kUnbounded;
+  const bool has_is = this->col_index_.size() || this->row_index_.size();
+  // If the model is known to be feasible, then there should be no IS,
+  // and all columns and rows are kIisStatusNotInConflict
+  if (is_feasible) assert(!has_is);
+  // If there is an IS, then all columns and rows not in the IS are
+  // kIisStatusNotInConflict
+  const HighsInt default_iis_status = is_feasible || has_is ?
+    kIisStatusNotInConflict :
+    kIisStatusMaybeInConflict;
+  return default_iis_status;
+}
+  
+void HighsIis::setStatus(const HighsLp& lp,
+			 const HighsModelStatus& model_status) {
   if (!this->valid_) return;
-  this->col_status_.assign(lp.num_col_, kIisStatusNotInConflict);
-  this->row_status_.assign(lp.num_row_, kIisStatusNotInConflict);
+  const HighsInt non_is_status = nonIsStatus(model_status);
+  this->col_status_.assign(lp.num_col_, non_is_status);
+  this->row_status_.assign(lp.num_row_, non_is_status);
   HighsInt iis_num_col = this->col_index_.size();
   HighsInt iis_num_row = this->row_index_.size();
   for (HighsInt iisCol = 0; iisCol < iis_num_col; iisCol++)
@@ -809,7 +827,8 @@ HighsStatus HighsIis::compute(const HighsLp& lp, const HighsOptions& options,
 
 bool indexStatusOkReturn(const bool return_value) { return return_value; }
 
-bool HighsIis::indexStatusOk(const HighsLp& lp) const {
+bool HighsIis::indexStatusOk(const HighsLp& lp,
+			     const HighsModelStatus& model_status) const {
   HighsInt num_col = lp.num_col_;
   HighsInt num_row = lp.num_row_;
   bool col_status_size_ok =
@@ -843,26 +862,30 @@ bool HighsIis::indexStatusOk(const HighsLp& lp) const {
   // other cols and rows are kIisStatusNotConflict
   std::vector<HighsInt> col_status = col_status_;
   std::vector<HighsInt> row_status = row_status_;
+  const HighsInt illegal_status = -99;
   for (HighsInt iX = 0; iX < num_iis_col; iX++) {
     HighsInt iCol = this->col_index_[iX];
     if (col_status_[iCol] != true_iis ? kIisStatusInConflict
                                       : kIisStatusMaybeInConflict)
       return indexStatusOkReturn(false);
-    col_status[iCol] = -1;
+    col_status[iCol] = illegal_status;
   }
   for (HighsInt iX = 0; iX < num_iis_row; iX++) {
     HighsInt iRow = this->row_index_[iX];
     if (row_status_[iRow] != true_iis ? kIisStatusInConflict
                                       : kIisStatusMaybeInConflict)
       return indexStatusOkReturn(false);
-    row_status[iRow] = -1;
+    row_status[iRow] = illegal_status;
   }
+  const HighsInt non_is_status = nonIsStatus(model_status);
   for (HighsInt iCol = 0; iCol < num_col; iCol++) {
-    if (col_status[iCol] >= 0 && col_status[iCol] != kIisStatusNotInConflict)
+    if (col_status[iCol] > illegal_status &&
+	col_status[iCol] != non_is_status)
       return indexStatusOkReturn(false);
   }
   for (HighsInt iRow = 0; iRow < num_row; iRow++) {
-    if (row_status[iRow] >= 0 && row_status[iRow] != kIisStatusNotInConflict)
+    if (row_status[iRow] > illegal_status &&
+	row_status[iRow] != non_is_status)
       return indexStatusOkReturn(false);
   }
   return indexStatusOkReturn(true);

highs/lp_data/HighsIis.h

@@ -44,15 +44,18 @@ class HighsIis {
   HighsStatus deduce(const HighsLp& lp, const HighsOptions& options,
                      const HighsBasis& basis);
   void setLp(const HighsLp& lp);
-  void setStatus(const HighsLp& lp);
+  HighsInt nonIsStatus(const HighsModelStatus& model_status) const;
+  void setStatus(const HighsLp& lp,
+		 const HighsModelStatus& model_status);
 
   HighsStatus compute(const HighsLp& lp, const HighsOptions& options,
                       const HighsBasis* basis = nullptr);
 
   bool trivial(const HighsLp& lp, const HighsOptions& options);
   bool rowValueBounds(const HighsLp& lp, const HighsOptions& options);
 
-  bool indexStatusOk(const HighsLp& lp) const;
+  bool indexStatusOk(const HighsLp& lp,
+		     const HighsModelStatus& model_status) const;
   bool lpDataOk(const HighsLp& lp, const HighsOptions& options) const;
   bool lpOk(const HighsOptions& options) const;
 

highs/lp_data/HighsInterface.cpp

@@ -1870,9 +1870,9 @@ HighsStatus Highs::getIisInterfaceReturn(const HighsStatus return_status) {
   }
   // Construct the ISS status vectors for cols and rows of original
   // model
-  this->iis_.setStatus(lp);
+  this->iis_.setStatus(lp, this->model_status_);
   // Check consistency of the col/row_index_ and col/row_status_
-  bool index_status_ok = this->iis_.indexStatusOk(lp);
+  bool index_status_ok = this->iis_.indexStatusOk(lp, this->model_status_);
   assert(index_status_ok);
   if (!index_status_ok) return HighsStatus::kError;
 
@@ -2102,14 +2102,14 @@ HighsStatus Highs::elasticityFilterReturn(
     getKktFailures(options_, model_, solution_, basis_, info_);
     info_.valid = true;
   }
-
-  // If the model is feasible, then the status of model is not known
-  if (feasible_model) this->model_status_ = HighsModelStatus::kNotset;
-  // The
-  assert(!feasible_model == (model_status == HighsModelStatus::kInfeasible));
   // The elasticity filter may well have identified infeasiblility, so
   // set this->model_status_
-  if (!feasible_model) this->model_status_ = HighsModelStatus::kInfeasible;
+  if (!feasible_model) {
+    this->model_status_ = HighsModelStatus::kInfeasible;
+  } else {
+    // If the model is feasible, then the status of model is not known
+    this->model_status_ = HighsModelStatus::kNotset;
+  }
   this->iis_ = iis;
   return return_status;
 }