Added Highs::getFixedLp(HighsLp& lp)

Author: jajhall

check/TestMipSolver.cpp

@@ -1035,3 +1035,41 @@ TEST_CASE("issue-2432", "[highs_test_mip_solver]") {
         "found\n");
   solve(highs, kHighsOffString, require_model_status, optimal_objective);
 }
+
+TEST_CASE("get-fixed-lp", "[highs_test_mip_solver]") {
+  std::string model = "avgas";
+  std::string model_file =
+      std::string(HIGHS_DIR) + "/check/instances/" + model + ".mps";
+  HighsLp fixed_lp;
+  Highs h;
+  h.setOptionValue("output_flag", dev_run);
+  REQUIRE(h.readModel(model_file) == HighsStatus::kOk);
+  REQUIRE(h.getFixedLp(fixed_lp) == HighsStatus::kError);
+
+  model = "flugpl";
+  model_file = std::string(HIGHS_DIR) + "/check/instances/" + model + ".mps";
+  REQUIRE(h.readModel(model_file) == HighsStatus::kOk);
+  REQUIRE(h.getFixedLp(fixed_lp) == HighsStatus::kError);
+
+  REQUIRE(h.run() == HighsStatus::kOk);
+  double mip_optimal_objective = h.getInfo().objective_function_value;
+  HighsSolution solution = h.getSolution();
+
+  REQUIRE(h.getFixedLp(fixed_lp) == HighsStatus::kOk);
+
+  REQUIRE(h.passModel(fixed_lp) == HighsStatus::kOk);
+  REQUIRE(h.run() == HighsStatus::kOk);
+
+  REQUIRE(h.getInfo().objective_function_value == mip_optimal_objective);
+
+  // Now run from saved solution without presolve
+  h.clearSolver();
+  h.setSolution(solution);
+  h.setOptionValue("presolve", kHighsOffString);
+  REQUIRE(h.run() == HighsStatus::kOk);
+
+  REQUIRE(h.getInfo().objective_function_value == mip_optimal_objective);
+  REQUIRE(h.getInfo().simplex_iteration_count == 0);
+
+  h.resetGlobalScheduler(true);
+}

highs/Highs.h

@@ -459,12 +459,11 @@ class Highs {
   }
 
   /**
-   * @brief Return a const pointer to the original row indices for the
-   * presolved model
+   * @brief Return an LP associated with a MIP and its solution, with
+   * each integer variable fixed to the value it takes in the MIP
+   * solution. If no solution is available, an error is returned.
    */
-  const HighsInt* getPresolveOrigRowsIndex() const {
-    return presolve_.data_.postSolveStack.getOrigRowsIndex();
-  }
+  HighsStatus getFixedLp(HighsLp& lp) const;
 
   /**
    * @brief Return a const reference to the incumbent LP

highs/lp_data/Highs.cpp

@@ -1893,6 +1893,47 @@ HighsStatus Highs::getStandardFormLp(HighsInt& num_col, HighsInt& num_row,
   return HighsStatus::kOk;
 }
 
+HighsStatus Highs::getFixedLp(HighsLp& lp) const {
+  if (!this->model_.lp_.isMip()) {
+    highsLogUser(options_.log_options, HighsLogType::kError,
+                 "Incumbent model is not a MIP, so cannot form fixed LP\n");
+    return HighsStatus::kError;
+  }
+  if (!this->solution_.value_valid) {
+    highsLogUser(options_.log_options, HighsLogType::kError,
+                 "Incumbent model does not have a valid solution, so cannot "
+                 "form fixed LP\n");
+    return HighsStatus::kError;
+  }
+  lp = this->model_.lp_;
+  const std::vector<HighsVarType> integrality = this->model_.lp_.integrality_;
+  lp.integrality_.clear();
+  HighsInt num_non_conts_fractional = 0;
+  double max_fractional = 0;
+  for (HighsInt iCol = 0; iCol < this->model_.lp_.num_col_; iCol++) {
+    if (integrality[iCol] == HighsVarType::kContinuous) continue;
+    double value = this->solution_.col_value[iCol];
+    if (integrality[iCol] == HighsVarType::kInteger ||
+        integrality[iCol] == HighsVarType::kSemiInteger) {
+      double fractional = fractionality(value);
+      if (fractional > this->options_.mip_feasibility_tolerance) {
+        num_non_conts_fractional++;
+        max_fractional = std::max(fractional, max_fractional);
+      }
+    }
+    lp.col_lower_[iCol] = value;
+    lp.col_upper_[iCol] = value;
+  }
+  if (num_non_conts_fractional) {
+    highsLogUser(
+        options_.log_options, HighsLogType::kWarning,
+        "Fixed LP has %d variables fixed at max fractional value of %g\n",
+        int(num_non_conts_fractional), max_fractional);
+    return HighsStatus::kWarning;
+  }
+  return HighsStatus::kOk;
+}
+
 HighsStatus Highs::getDualRay(bool& has_dual_ray, double* dual_ray_value) {
   has_dual_ray = false;
   return getDualRayInterface(has_dual_ray, dual_ray_value);