Testing Ines, and now calling FJ after trivial heuristics and FJ

Author: jajhall

highs/lp_data/HighsOptions.h

@@ -439,6 +439,7 @@ struct HighsOptionsStruct {
   double mip_rel_gap;
   double mip_abs_gap;
   double mip_heuristic_effort;
+  bool mip_heuristic_run_ines;
   bool mip_heuristic_run_feasibility_jump;
   bool mip_heuristic_run_rins;
   bool mip_heuristic_run_rens;
@@ -590,6 +591,7 @@ struct HighsOptionsStruct {
         mip_rel_gap(0.0),
         mip_abs_gap(0.0),
         mip_heuristic_effort(0.0),
+        mip_heuristic_run_ines(false),
         mip_heuristic_run_feasibility_jump(false),
         mip_heuristic_run_rins(false),
         mip_heuristic_run_rens(false),
@@ -1136,6 +1138,11 @@ class HighsOptions : public HighsOptionsStruct {
         &mip_heuristic_effort, 0.0, 0.05, 1.0);
     records.push_back(record_double);
 
+    record_bool =
+        new OptionRecordBool("mip_heuristic_run_ines", "Use Ines' heuristic",
+                             advanced, &mip_heuristic_run_ines, true);
+    records.push_back(record_bool);
+
     record_bool =
         new OptionRecordBool("mip_heuristic_run_feasibility_jump",
                              "Use the feasibility jump heuristic", advanced,

highs/mip/HighsFeasibilityJump.cpp

@@ -10,6 +10,7 @@
 #include "util/HighsSparseMatrix.h"
 
 HighsModelStatus HighsMipSolverData::feasibilityJump() {
+  this->feasibility_jump_time_ = -mipsolver.timer_.read();
   // This is the (presolved) model being solved
   const HighsLp* model = this->mipsolver.model_;
   const HighsLogOptions& log_options = mipsolver.options_mip_->log_options;
@@ -35,6 +36,9 @@ HighsModelStatus HighsMipSolverData::feasibilityJump() {
       /* equalityTolerance = */ epsilon,
       /* violationTolerance = */ feastol);
 
+  bool have_solution = this->mipsolver.solution_.size() > 0;
+  if (have_solution) col_value = this->mipsolver.solution_;
+
   for (HighsInt col = 0; col < model->num_col_; ++col) {
     double lower = model->col_lower_[col];
     double upper = model->col_upper_[col];
@@ -66,13 +70,23 @@ HighsModelStatus HighsMipSolverData::feasibilityJump() {
     solver.addVar(fjVarType, lower, upper,
                   sense_multiplier * model->col_cost_[col]);
 
-    double initial_assignment = 0.0;
-    if (std::isfinite(lower)) {
-      initial_assignment = lower;
-    } else if (std::isfinite(upper)) {
-      initial_assignment = upper;
+    if (have_solution) {
+      if (col_value[col] < lower || col_value[col] > upper) {
+        highsLogUser(log_options, HighsLogType::kError,
+                     "HighsMipSolverData::feasibilityJump() incumbent solution "
+                     "for column %d is %g, not in [%g, %g]\n",
+                     int(col), col_value[col], lower, upper);
+        return HighsModelStatus::kSolveError;
+      }
+    } else {
+      double initial_assignment = 0.0;
+      if (std::isfinite(lower)) {
+        initial_assignment = lower;
+      } else if (std::isfinite(upper)) {
+        initial_assignment = upper;
+      }
+      col_value[col] = initial_assignment;
     }
-    col_value[col] = initial_assignment;
   }
 
   HighsSparseMatrix a_matrix;
@@ -129,6 +143,8 @@ HighsModelStatus HighsMipSolverData::feasibilityJump() {
     // use trySolution rather than addIncumbent for an explicit check.
     trySolution(col_value, kSolutionSourceFeasibilityJump);
   }
+  this->feasibility_jump_time_ += mipsolver.timer_.read();
+  this->feasibility_jump_objective_ = objective_function_value;
   return HighsModelStatus::kNotset;
 #endif
 }

highs/mip/HighsMipSolver.cpp

@@ -183,11 +183,31 @@ void HighsMipSolver::run() {
       mipdata_->callbackUserSolution(solution_objective_,
                                      kUserMipSolutionCallbackOriginAfterSetup);
 
-    const bool mip_heuristic_run_ines = false;
+    // Apply the trivial heuristics
+    analysis_.mipTimerStart(kMipClockTrivialHeuristics);
+    HighsModelStatus returned_model_status = mipdata_->trivialHeuristics();
+    analysis_.mipTimerStop(kMipClockTrivialHeuristics);
+    if (modelstatus_ == HighsModelStatus::kNotset &&
+        returned_model_status == HighsModelStatus::kInfeasible) {
+      // trivialHeuristics can spot trivial infeasibility, so act on it
+      modelstatus_ = returned_model_status;
+      cleanupSolve();
+      return;
+    }
+
+    // Don't report on this MIP type, or re-record the MIP
     const bool ines_silent = true;
     if (mipdata_->mipIsInes(ines_silent)) {
-      if (mip_heuristic_run_ines) {
+      assert(mipdata_->mipIsInesLookup());
+      if (options_mip_->mip_heuristic_run_ines) {
+        analysis_.mipTimerStart(kMipClockInes);
         HighsStatus ines_status = mipdata_->heuristics.mipHeuristicInes();
+        analysis_.mipTimerStop(kMipClockInes);
+        if (ines_status == HighsStatus::kError) {
+          modelstatus_ = HighsModelStatus::kSolveError;
+          cleanupSolve();
+          return;
+        }
       }
     }
 
@@ -203,25 +223,27 @@ void HighsMipSolver::run() {
         cleanupSolve();
         return;
       }
-      const bool bailout_after_feasibility_jump = false;
-      if (bailout_after_feasibility_jump) {
-        highsLogUser(options_mip_->log_options, HighsLogType::kInfo,
-                     "HighsMipSolver: Bailing out after Feasibility Jump with "
-                     "model status = %s\n",
-                     utilModelStatusToString(returned_model_status).c_str());
-        modelstatus_ = HighsModelStatus::kInterrupt;
+      if (returned_model_status == HighsModelStatus::kSolveError) {
+        modelstatus_ = returned_model_status;
         cleanupSolve();
         return;
       }
     }
-    // Apply the trivial heuristics
-    analysis_.mipTimerStart(kMipClockTrivialHeuristics);
-    HighsModelStatus returned_model_status = mipdata_->trivialHeuristics();
-    analysis_.mipTimerStop(kMipClockTrivialHeuristics);
-    if (modelstatus_ == HighsModelStatus::kNotset &&
-        returned_model_status == HighsModelStatus::kInfeasible) {
-      // trivialHeuristics can spot trivial infeasibility, so act on it
-      modelstatus_ = returned_model_status;
+    if (!submip && mipdata_->mipIsInesLookup()) {
+      printf("grepInesFJ,%s,%d,%d,%g,%d,%g,%g,%g,%g\n",
+             model_->model_name_.c_str(), int(model_->num_col_),
+             int(model_->num_row_), model_->offset_,
+             mipdata_->ines_zero_solution_, mipdata_->ines_objective_,
+             mipdata_->ines_time_, mipdata_->feasibility_jump_objective_,
+             mipdata_->feasibility_jump_time_);
+    }
+    const bool bailout_after_feasibility_jump = false;
+    if (bailout_after_feasibility_jump) {
+      highsLogUser(options_mip_->log_options, HighsLogType::kInfo,
+                   "HighsMipSolver: Bailing out after Feasibility Jump with "
+                   "model status = %s\n",
+                   utilModelStatusToString(returned_model_status).c_str());
+      modelstatus_ = HighsModelStatus::kInterrupt;
       cleanupSolve();
       return;
     }

highs/mip/HighsMipSolverData.cpp

@@ -81,6 +81,9 @@ std::string HighsMipSolverData::solutionSourceToString(
   } else if (solution_source == kSolutionSourceUserSolution) {
     if (code) return "X";
     return "User solution";
+  } else if (solution_source == kSolutionSourceInes) {
+    if (code) return "i";
+    return "Ines";
   } else if (solution_source == kSolutionSourceCleanup) {
     if (code) return " ";
     return "";
@@ -687,6 +690,11 @@ void HighsMipSolverData::init() {
   knapsack_capacity_ = 0;
   knapsack_integral_scale_ = 0;
   mip_problem_data_.clear();
+  ines_time_ = 0;
+  ines_objective_ = 0;
+  ines_zero_solution_ = false;
+  feasibility_jump_time_ = 0;
+  feasibility_jump_objective_ = 0;
 
   if (mipsolver.options_mip_->mip_report_level == 0)
     dispfreq = 0;
@@ -2730,6 +2738,13 @@ bool HighsMipSolverData::mipIsInes(const bool silent) {
   return true;
 }
 
+bool HighsMipSolverData::mipIsInesLookup() const {
+  HighsInt num_mip = static_cast<HighsInt>(mip_problem_data_.record.size());
+  if (num_mip == 0) return false;
+  return mip_problem_data_.record[num_mip - 1].type ==
+         HighsMipProblemType::kInes;
+}
+
 void HighsMipSolverData::mipIsOther() {
   const HighsLp& lp = *(mipsolver.model_);
   HighsMipProblemRecord mip_problem_record;

highs/mip/HighsMipSolverData.h

@@ -85,6 +85,7 @@ enum MipSolutionSource : int {
   kSolutionSourceTrivialL,
   kSolutionSourceTrivialU,
   kSolutionSourceTrivialP,
+  kSolutionSourceInes,
   kSolutionSourceCleanup,
   kSolutionSourceCount
 };
@@ -181,6 +182,11 @@ struct HighsMipSolverData {
   HighsInt knapsack_integral_scale_;
 
   HighsMipProblemData mip_problem_data_;
+  double ines_time_;
+  double ines_objective_;
+  bool ines_zero_solution_;
+  double feasibility_jump_time_;
+  double feasibility_jump_objective_;
 
   HighsMipSolverData(HighsMipSolver& mipsolver)
       : mipsolver(mipsolver),
@@ -331,6 +337,7 @@ struct HighsMipSolverData {
       const double mipsolver_objective_value,
       const userMipSolutionCallbackOrigin user_solution_callback_origin);
   bool mipIsKnapsack(const bool silent = false);
+  bool mipIsInesLookup() const;
   bool mipIsInes(const bool silent = false);
   void mipIsOther();
 };