Don't pass test2 of AMPL infeasible MIP

Author: jajhall

check/TestIis.cpp

@@ -472,6 +472,8 @@ TEST_CASE("lp-get-iis-avgas", "[iis]") {
 
 TEST_CASE("lp-feasibility-relaxation", "[iis]") {
   // Using infeasible LP from AMPL documentation
+  //
+  // https://mp.ampl.com/features-guide.html#feasibility-relaxation
   HighsLp lp;
   lp.model_name_ = "ampl_infeas";
   lp.num_col_ = 2;
@@ -488,14 +490,14 @@ TEST_CASE("lp-feasibility-relaxation", "[iis]") {
   lp.a_matrix_.value_ = {-1, -3, 20, 21, 2, 1};
   lp.integrality_ = {HighsVarType::kInteger, HighsVarType::kInteger};
   Highs h;
-  h.setOptionValue("output_flag", dev_run);
+  //h.setOptionValue("output_flag", dev_run);
   const HighsSolution& solution = h.getSolution();
   h.passModel(lp);
 
-  const bool all_tests = false;
+  const bool all_tests = false;//true;
   const bool test0 = false || all_tests;
   const bool test1 = false || all_tests;
-  const bool test2 = false || all_tests;
+  const bool test2 = true || all_tests;
   const bool test3 = false || all_tests;
   if (test0) {
     // Vanilla feasibility relaxation
@@ -531,7 +533,8 @@ TEST_CASE("lp-feasibility-relaxation", "[iis]") {
           "===============================\n"
           "Local penalties RHS {1, -1, 10}\n"
           "===============================\n");
-    // Now test local penalties
+    // Now test local penalties, allowing lower bounds to be violated,
+    // but not upper bounds
     //
     // constraint 0: normal weight
     //
@@ -544,7 +547,7 @@ TEST_CASE("lp-feasibility-relaxation", "[iis]") {
                             local_rhs_penalty.data());
     // Should get slacks (-3, 4, 0)
     h.writeSolution("", 1);
-    REQUIRE(solution.row_value[0] == lp.row_lower_[0] - 3);
+    REQUIRE(solution.row_value[0] == lp.row_lower_[0] + 3);
     REQUIRE(solution.row_value[1] == lp.row_upper_[1] - 4);
     REQUIRE(solution.row_value[2] == lp.row_upper_[2]);
   }
@@ -566,6 +569,7 @@ TEST_CASE("lp-feasibility-relaxation", "[iis]") {
     h.feasibilityRelaxation(1, 1, 0, nullptr, nullptr,
                             local_rhs_penalty.data());
     // Should get slacks (18, 2, -1)
+    h.writeSolution("", 1);
     REQUIRE(solution.row_value[0] == lp.row_lower_[0] + 18);
     REQUIRE(solution.row_value[1] == lp.row_upper_[1] - 2);
     REQUIRE(solution.row_value[2] == lp.row_upper_[2] + 1);

highs/Highs.h

@@ -1728,7 +1728,7 @@ class Highs {
   HighsStatus getIisInterfaceReturn(const HighsStatus return_status);
 
   HighsStatus elasticityFilterReturn(
-      const HighsStatus return_status, const bool feasible_model,
+      const HighsStatus return_status, 
       const std::string& original_model_name, const HighsInt original_num_col,
       const HighsInt original_num_row,
       const std::vector<double>& original_col_cost,

highs/lp_data/HighsInterface.cpp

@@ -2047,7 +2047,7 @@ HighsStatus Highs::getIisInterface() {
 }
 
 HighsStatus Highs::elasticityFilterReturn(
-    const HighsStatus return_status, const bool feasible_model,
+    const HighsStatus return_status, 
     const std::string& original_model_name, const HighsInt original_num_col,
     const HighsInt original_num_row,
     const std::vector<double>& original_col_cost,
@@ -2059,7 +2059,11 @@ HighsStatus Highs::elasticityFilterReturn(
   // LP, so save them
   const bool iis_valid = this->iis_.valid_;
   HighsIis iis;
-  if (iis_valid) iis = this->iis_;
+  if (iis_valid) {
+    iis = this->iis_;
+  } else {
+    assert(return_status != HighsStatus::kOk);
+  }
   double objective_function_value = info_.objective_function_value;
   // Delete any additional rows and columns, and restore the original
   // column costs and bounds
@@ -2108,10 +2112,17 @@ 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;
   // If there was valid HighsIis data, then restore it
-  if (iis_valid) this->iis_ = iis;
+  if (iis_valid) {
+    // If the model is feasible, then the status of model is not known
+    if (iis.status_ == kIisModelStatusFeasible)
+      this->model_status_ = HighsModelStatus::kNotset;
+    this->iis_ = iis;
+  } else {
+    // No valid HighsIis should imply a runtime error, so ensure that
+    // no misleading model status is set
+    this->model_status_ = HighsModelStatus::kNotset;
+  }
   return return_status;
 }
 
@@ -2434,20 +2445,24 @@ HighsStatus Highs::elasticityFilter(const double global_lower_penalty,
     return run_status;
   };
 
+  // Solve the elastic LP
   run_status = solveLp();
 
   if (run_status != HighsStatus::kOk)
-    return elasticityFilterReturn(run_status, false, original_model_name,
+    return elasticityFilterReturn(run_status, original_model_name,
                                   original_num_col, original_num_row,
                                   original_col_cost, original_col_lower,
                                   original_col_upper, original_integrality);
   if (kIisDevReport) this->writeSolution("", kSolutionStylePretty);
   // Model status should be optimal, unless model is unbounded
   assert(this->model_status_ == HighsModelStatus::kOptimal ||
          this->model_status_ == HighsModelStatus::kUnbounded);
-
-  if (!get_iis)
-    return elasticityFilterReturn(HighsStatus::kOk, false, original_model_name,
+  this->iis_.valid_ = true;
+  this->iis_.status_ = this->info_.objective_function_value > 0 ?
+    kIisModelStatusInfeasible :
+    kIisModelStatusFeasible;
+  if (!get_iis) 
+    return elasticityFilterReturn(HighsStatus::kOk, original_model_name,
                                   original_num_col, original_num_row,
                                   original_col_cost, original_col_lower,
                                   original_col_upper, original_integrality);
@@ -2458,7 +2473,6 @@ HighsStatus Highs::elasticityFilter(const double global_lower_penalty,
   const HighsSolution& solution = this->getSolution();
   // Now fix e-variables that are positive and re-solve until e-LP is infeasible
   HighsInt loop_k = 0;
-  bool feasible_model = false;
   for (;;) {
     if (kIisDevReport)
       printf("\nElasticity filter pass %d\n==============\n", int(loop_k));
@@ -2501,13 +2515,13 @@ HighsStatus Highs::elasticityFilter(const double global_lower_penalty,
     }
     if (num_fixed == 0) {
       // No elastic variables were positive, so problem is feasible
-      feasible_model = true;
+      this->iis_.status_ = kIisModelStatusFeasible;
       break;
     }
     HighsStatus run_status = solveLp();
     if (run_status != HighsStatus::kOk)
       return elasticityFilterReturn(
-          run_status, feasible_model, original_model_name, original_num_col,
+          run_status, original_model_name, original_num_col,
           original_num_row, original_col_cost, original_col_lower,
           original_col_upper, original_integrality);
     if (kIisDevReport) this->writeSolution("", kSolutionStylePretty);
@@ -2553,7 +2567,7 @@ HighsStatus Highs::elasticityFilter(const double global_lower_penalty,
     }
   }
   HighsInt num_iis_row = iis.row_index_.size();
-  if (feasible_model) {
+  if (iis.status_ == kIisModelStatusFeasible) {
     assert(num_enforced_col_ecol == 0 && num_enforced_row_ecol == 0);
     assert(num_iis_row == 0);
   }
@@ -2571,17 +2585,15 @@ HighsStatus Highs::elasticityFilter(const double global_lower_penalty,
 
   iis.valid_ = true;
   iis.strategy_ = this->options_.iis_strategy;
-  if (feasible_model) {
-    iis.status_ = kIisModelStatusFeasible;
+  if (iis.status_ == kIisModelStatusFeasible) 
     return elasticityFilterReturn(
-        HighsStatus::kOk, feasible_model, original_model_name, original_num_col,
+        HighsStatus::kOk, original_model_name, original_num_col,
         original_num_row, original_col_cost, original_col_lower,
         original_col_upper, original_integrality);
-  }
-  assert(!feasible_model);
   // Model is infeasible because there are (at least) a positive
   // number of rows in the infeasibility set. Hence the IIS status is
   // reducible
+  assert(this->iis_.status_ > kIisModelStatusFeasible);
   assert(num_iis_row > 0);
   iis.status_ = kIisModelStatusReducible;
 
@@ -2652,7 +2664,7 @@ HighsStatus Highs::elasticityFilter(const double global_lower_penalty,
   assert(iis.row_bound_.size() == iis.row_index_.size());
   for (HighsInt iX = 0; iX < num_iis_row; iX++) assert(iis.row_bound_[iX] >= 0);
   return elasticityFilterReturn(
-      HighsStatus::kOk, feasible_model, original_model_name, original_num_col,
+      HighsStatus::kOk, original_model_name, original_num_col,
       original_num_row, original_col_cost, original_col_lower,
       original_col_upper, original_integrality);
 }