Made @mathgeekcoder's changes, and computing the correct reduced objective value for maximization, as well as minimization problems; formatted

Author: jajhall

check/TestMipSolver.cpp

@@ -1037,8 +1037,8 @@ TEST_CASE("issue-2432", "[highs_test_mip_solver]") {
 }
 
 TEST_CASE("mip-race", "[highs_test_mip_solver]") {
-  const bool ci_test = false;
-  const std::string test_build_model = "fiball";
+  const bool ci_test = true;
+  const std::string test_build_model = "neos-3381206-awhea";  //"fiball";
   const std::string model = ci_test ? "flugpl" : test_build_model;
   // "neos-3381206-awhea";
   const std::string model_file =
@@ -1050,8 +1050,18 @@ TEST_CASE("mip-race", "[highs_test_mip_solver]") {
   h.setOptionValue("mip_race_concurrency", mip_race_concurrency);
   h.setOptionValue("mip_race_read_solutions", true);
   REQUIRE(h.readModel(model_file) == HighsStatus::kOk);
-  REQUIRE(h.run() == HighsStatus::kOk);
-  REQUIRE(h.getModelStatus() == HighsModelStatus::kOptimal);
+  for (int k = 0; k < 2; k++) {
+    if (k == 1) {
+      HighsLp lp = h.getLp();
+      for (HighsInt iCol = 0; iCol < lp.num_col_; iCol++)
+        lp.col_cost_[iCol] = -lp.col_cost_[iCol];
+      REQUIRE(h.changeColsCost(0, lp.num_col_ - 1, lp.col_cost_.data()) ==
+              HighsStatus::kOk);
+      REQUIRE(h.changeObjectiveSense(ObjSense::kMaximize) == HighsStatus::kOk);
+    }
+    REQUIRE(h.run() == HighsStatus::kOk);
+    REQUIRE(h.getModelStatus() == HighsModelStatus::kOptimal);
+  }
 
   if (ci_test) {
     h.clearSolver();

highs/Highs.h

@@ -1717,7 +1717,7 @@ class Highs {
   HighsStatus mipRaceResults(HighsMipSolverInfo& mip_solver_info,
                              const std::vector<HighsMipSolverInfo>& worker_info,
                              const std::vector<double>& mip_time,
-			     const double& report_mip_time);
+                             const double& report_mip_time);
 };
 
 // Start of deprecated methods not in the Highs class

highs/lp_data/Highs.cpp

@@ -4106,9 +4106,8 @@ HighsStatus Highs::callSolveMip() {
         });
     loop_mip_time += timer_.read();
     // Determine the winner and report on the solution
-    HighsStatus call_status =
-      this->mipRaceResults(mip_solver_info, worker_info, mip_time,
-			   loop_mip_time);
+    HighsStatus call_status = this->mipRaceResults(mip_solver_info, worker_info,
+                                                   mip_time, loop_mip_time);
     if (call_status == HighsStatus::kError) {
       const bool undo_mods = true;
       return returnFromOptimizeModel(HighsStatus::kError, undo_mods);

highs/lp_data/HighsInterface.cpp

@@ -4264,8 +4264,7 @@ void HighsMipSolverInfo::clear() {
 HighsStatus Highs::mipRaceResults(
     HighsMipSolverInfo& mip_solver_info,
     const std::vector<HighsMipSolverInfo>& worker_info,
-    const std::vector<double>& mip_time,
-    const double& report_mip_time) {
+    const std::vector<double>& mip_time, const double& report_mip_time) {
   const HighsInt mip_race_concurrency = this->options_.mip_race_concurrency;
   HighsInt winning_instance = -1;
   HighsModelStatus winning_model_status = HighsModelStatus::kNotset;
@@ -4276,7 +4275,7 @@ HighsStatus Highs::mipRaceResults(
         instance == 0 ? mip_solver_info : worker_info[instance];
     HighsModelStatus instance_model_status = solver_info.modelstatus;
     highsLogUser(options_.log_options, HighsLogType::kInfo,
-                 "   Solver %2d has best objective %15.8g, gap %6.2f\% (time "
+                 "   Solver %2d has best objective %15.8g, gap %6.2f%% (time "
                  "= %6.2f), and status %s\n",
                  int(instance), solver_info.solution_objective,
                  1e2 * solver_info.gap, mip_time[instance],

highs/mip/HighsMipSolver.h

@@ -8,6 +8,8 @@
 #ifndef MIP_HIGHS_MIP_SOLVER_H_
 #define MIP_HIGHS_MIP_SOLVER_H_
 
+#include <atomic>
+
 #include "Highs.h"
 #include "lp_data/HighsCallback.h"
 #include "lp_data/HighsOptions.h"
@@ -21,6 +23,7 @@ class HighsImplications;
 
 const HighsInt kMipRaceNoSolution = -1;
 
+/*
 struct MipRaceIncumbent {
   HighsInt start_write_incumbent = kMipRaceNoSolution;
   HighsInt finish_write_incumbent = kMipRaceNoSolution;
@@ -32,11 +35,11 @@ struct MipRaceIncumbent {
   HighsInt read(const HighsInt last_incumbent_read, double& objective_,
                 std::vector<double>& solution_) const;
 };
+*/
 
-/*
-  struct MipRaceIncumbent {
-  std::atomic<HighsInt> start_write_incumbent = kMipRaceNoSolution;
-  std::atomic<HighsInt> finish_write_incumbent = kMipRaceNoSolution;
+struct MipRaceIncumbent {
+  std::atomic<HighsInt> start_write_incumbent{kMipRaceNoSolution};
+  std::atomic<HighsInt> finish_write_incumbent{kMipRaceNoSolution};
   double objective = -kHighsInf;
   std::vector<double> solution;
   void clear();
@@ -61,7 +64,7 @@ struct MipRaceIncumbent {
     solution = std::move(moving.solution);
   }
 };
-*/
+
 struct MipRaceRecord {
   std::vector<MipRaceIncumbent> incumbent;
   void clear();
@@ -83,7 +86,7 @@ struct MipRace {
                   const HighsLogOptions log_options_);
   HighsInt concurrency() const;
   void update(const double objective, const std::vector<double>& solution);
-  bool newSolution(const HighsInt instance, double objective,
+  bool newSolution(const HighsInt instance, double& objective,
                    std::vector<double>& solution);
   void report() const;
 };

highs/mip/HighsMipSolverData.cpp

@@ -1401,6 +1401,20 @@ bool HighsMipSolverData::addIncumbent(const std::vector<double>& sol,
                                      sol, possibly_store_as_new_incumbent)
                                : 0;
 
+  if (solution_source == kSolutionSourceHighsSolution) {
+    printf(
+        "HighsMipSolverData::addIncumbent HiGHS solution Offset = %15.8g; Obj "
+        "= %15.8g; UB = %15.8g; PossAdd = %s",
+        mipsolver.model_->offset_, solobj, upper_bound,
+        possibly_store_as_new_incumbent ? "T" : "F");
+    if (possibly_store_as_new_incumbent) {
+      printf("; TransObj = %15.8g; TransSolobj < UB %s \n", transformed_solobj,
+             transformed_solobj < upper_bound ? "T" : "F");
+    } else {
+      printf("\n");
+    }
+    fflush(stdout);
+  }
   if (possibly_store_as_new_incumbent) {
     solobj = transformed_solobj;
     if (solobj >= upper_bound) return false;
@@ -2677,18 +2691,57 @@ void HighsMipSolverData::queryExternalSolution(
   if (!mipsolver.options_mip_->mip_race_read_solutions) return;
   MipRace& mip_race = mipsolver.mip_race_;
   if (!mip_race.record) return;
-  double instance_solution_objective_value = kHighsInf;
+  double instance_objective_value = kHighsInf;
   std::vector<double> instance_solution;
   for (HighsInt instance = 0; instance < mipRaceConcurrency(); instance++) {
     if (instance == mip_race.my_instance) continue;
-    if (!mip_race.newSolution(instance, instance_solution_objective_value,
+    if (!mip_race.newSolution(instance, instance_objective_value,
                               instance_solution))
       continue;
     // Have read a new incumbent
-    std::vector<double> reduced_instance_solution;
-    reduced_instance_solution =
+    //
+    // Objective is assumed to be original_offset + (original_c)^T(original_x),
+    // but MIP solver bounds are based on the reduced objective
+    // (reduced_c)^T(reduced_x)
+    //
+    // Now, original_sense*[reduced_offset + (reduced_c)^T(reduced_x)] is an
+    // objective in the original space, so
+    //
+    // f0 + c0^Tx0 = s*(f1 + c1^Tx1)
+    //
+    // where 0 => original; 1 => reduced
+    //
+    // This allows the reduced objective value to be deduced as
+    //
+    // c1^Tx1 = s*(f0 + c0^Tx0) - f1
+    //
+    // (reduced_c)^T(reduced_x) = original_sense*[original_offset +
+    // (original_c)^T(original_x) - reduced_offset]
+    //
+    double reduced_instance_objective_value = instance_objective_value;
+    reduced_instance_objective_value *= int(mipsolver.orig_model_->sense_);
+    reduced_instance_objective_value -= mipsolver.model_->offset_;
+    // Get the solution in the reduced space
+    std::vector<double> reduced_instance_solution =
         postSolveStack.getReducedPrimalSolution(instance_solution);
-    addIncumbent(reduced_instance_solution, instance_solution_objective_value,
+
+    double check_objective_value = 0;
+    for (HighsInt iCol = 0; iCol < mipsolver.model_->num_col_; iCol++)
+      check_objective_value +=
+          mipsolver.colCost(iCol) * reduced_instance_solution[iCol];
+    double dl_objective_value =
+        std::fabs(check_objective_value - reduced_instance_objective_value);
+    assert(dl_objective_value < 1e-12 * (1 + std::fabs(check_objective_value)));
+    printf(
+        "HighsMipSolverData::queryExternalSolution: (sense = %d; offset = "
+        "%11.4g) modified objective from %11.4g to %11.4g (Check = %11.4g; "
+        "Delta = %11.4g)\n",
+        int(mipsolver.orig_model_->sense_), mipsolver.model_->offset_,
+        instance_objective_value, reduced_instance_objective_value,
+        check_objective_value, dl_objective_value);
+    fflush(stdout);
+
+    addIncumbent(reduced_instance_solution, reduced_instance_objective_value,
                  kSolutionSourceHighsSolution);
   }
 }
@@ -2956,15 +3009,15 @@ void MipRaceRecord::report(const HighsLogOptions log_options) const {
   highsLogUser(log_options, HighsLogType::kInfo, "\nStartWrite:        ");
   for (HighsInt instance = 0; instance < mip_race_concurrency; instance++)
     highsLogUser(log_options, HighsLogType::kInfo, " %20d",
-                 this->incumbent[instance].start_write_incumbent);
+                 int(this->incumbent[instance].start_write_incumbent));
   highsLogUser(log_options, HighsLogType::kInfo, "\nObjective:         ");
   for (HighsInt instance = 0; instance < mip_race_concurrency; instance++)
     highsLogUser(log_options, HighsLogType::kInfo, " %20.12g",
                  this->incumbent[instance].objective);
   highsLogUser(log_options, HighsLogType::kInfo, "\nFinishWrite:       ");
   for (HighsInt instance = 0; instance < mip_race_concurrency; instance++)
     highsLogUser(log_options, HighsLogType::kInfo, " %20d",
-                 this->incumbent[instance].finish_write_incumbent);
+                 int(this->incumbent[instance].finish_write_incumbent));
   highsLogUser(log_options, HighsLogType::kInfo, "\n");
 }
 
@@ -2997,7 +3050,7 @@ void MipRace::update(const double objective,
   //  this->report();
 }
 
-bool MipRace::newSolution(const HighsInt instance, double objective,
+bool MipRace::newSolution(const HighsInt instance, double& objective,
                           std::vector<double>& solution) {
   assert(this->record);
   HighsInt new_incumbent_read = this->record->incumbent[instance].read(