Heuristic Improvements: balance between generation and improvement heuristics

benchmarks/linear_programming/cuopt/mip_test_instances.hpp

@@ -17,175 +17,11 @@
 #pragma once
 #include <string>
 #include <vector>
-std::vector<std::string> instances = {"30n20b8.mps",
-                                      "50v-10.mps",
-                                      "CMS750_4.mps",
-                                      "academictimetablesmall.mps",
-                                      "air05.mps",
-                                      "app1-1.mps",
-                                      "app1-2.mps",
-                                      "assign1-5-8.mps",
-                                      "atlanta-ip.mps",
-                                      "bab2.mps",
-                                      "bab6.mps",
-                                      "beasleyC3.mps",
-                                      "binkar10_1.mps",
-                                      "blp-ar98.mps",
-                                      "blp-ic98.mps",
-                                      "bppc4-08.mps",
-                                      "brazil3.mps",
-                                      "cmflsp50-24-8-8.mps",
-                                      "co-100.mps",
-                                      "cod105.mps",
-                                      "comp07-2idx.mps",
-                                      "comp21-2idx.mps",
-                                      "csched007.mps",
-                                      "csched008.mps",
-                                      "cvs16r128-89.mps",
-                                      "dano3_3.mps",
-                                      "decomp2.mps",
-                                      "drayage-100-23.mps",
-                                      "drayage-25-23.mps",
-                                      "eil33-2.mps",
-                                      "eilA101-2.mps",
-                                      "exp-1-500-5-5.mps",
-                                      "fast0507.mps",
-                                      "fastxgemm-n2r6s0t2.mps",
-                                      "fiball.mps",
-                                      "gen-ip002.mps",
-                                      "germanrr.mps",
-                                      "glass4.mps",
-                                      "graph20-20-1rand.mps",
-                                      "graphdraw-domain.mps",
-                                      "h80x6320d.mps",
-                                      "highschool1-aigio.mps",
-                                      "hypothyroid-k1.mps",
-                                      "icir97_tension.mps",
-                                      "irish-electricity.mps",
-                                      "istanbul-no-cutoff.mps",
-                                      "k1mushroom.mps",
-                                      "lectsched-5-obj.mps",
-                                      "leo1.mps",
-                                      "leo2.mps",
-                                      "lotsize.mps",
-                                      "mad.mps",
-                                      "map10.mps",
-                                      "map16715-04.mps",
-                                      "markshare2.mps",
-                                      "markshare_4_0.mps",
-                                      "mas74.mps",
-                                      "mc11.mps",
-                                      "mcsched.mps",
-                                      "mik-250-20-75-4.mps",
-                                      "momentum1.mps",
-                                      "mushroom-best.mps",
-                                      "mzzv11.mps",
-                                      "mzzv42z.mps",
-                                      "n2seq36q.mps",
-                                      "n3div36.mps",
-                                      "neos-1171448.mps",
-                                      "neos-1171737.mps",
-                                      "neos-1354092.mps",
-                                      "neos-1445765.mps",
-                                      "neos-1456979.mps",
-                                      "neos-1582420.mps",
-                                      "neos-2657525-crna.mps",
-                                      "neos-2746589-doon.mps",
-                                      "neos-3024952-loue.mps",
-                                      "neos-3046615-murg.mps",
-                                      "neos-3216931-puriri.mps",
-                                      "neos-3402294-bobin.mps",
-                                      "neos-3656078-kumeu.mps",
-                                      "neos-3754480-nidda.mps",
-                                      "neos-4300652-rahue.mps",
-                                      "neos-4338804-snowy.mps",
-                                      "neos-4387871-tavua.mps",
-                                      "neos-4413714-turia.mps",
-                                      "neos-4532248-waihi.mps",
-                                      "neos-4722843-widden.mps",
-                                      "neos-4738912-atrato.mps",
-                                      "neos-4763324-toguru.mps",
-                                      "neos-4954672-berkel.mps",
-                                      "neos-5049753-cuanza.mps",
-                                      "neos-5093327-huahum.mps",
-                                      "neos-5107597-kakapo.mps",
-                                      "neos-5114902-kasavu.mps",
-                                      "neos-5188808-nattai.mps",
-                                      "neos-5195221-niemur.mps",
-                                      "neos-662469.mps",
-                                      "neos-787933.mps",
-                                      "neos-848589.mps",
-                                      "neos-860300.mps",
-                                      "neos-911970.mps",
-                                      "neos-933966.mps",
-                                      "neos-950242.mps",
-                                      "neos17.mps",
-                                      "neos5.mps",
-                                      "net12.mps",
-                                      "netdiversion.mps",
-                                      "nexp-150-20-8-5.mps",
-                                      "ns1644855.mps",
-                                      "ns1760995.mps",
-                                      "ns1830653.mps",
-                                      "nursesched-medium-hint03.mps",
-                                      "nursesched-sprint02.mps",
-                                      "opm2-z10-s4.mps",
-                                      "pg.mps",
-                                      "physiciansched3-3.mps",
-                                      "piperout-08.mps",
-                                      "piperout-27.mps",
-                                      "pk1.mps",
-                                      "qap10.mps",
-                                      "radiationm18-12-05.mps",
-                                      "radiationm40-10-02.mps",
-                                      "rail01.mps",
-                                      "rail02.mps",
-                                      "rail507.mps",
-                                      "ran14x18-disj-8.mps",
-                                      "rmatr100-p10.mps",
-                                      "rmatr200-p5.mps",
-                                      "rocI-4-11.mps",
-                                      "rocII-5-11.mps",
-                                      "rococoB10-011000.mps",
-                                      "rococoC10-001000.mps",
-                                      "roi2alpha3n4.mps",
-                                      "roi5alpha10n8.mps",
-                                      "roll3000.mps",
-                                      "s100.mps",
-                                      "s250r10.mps",
-                                      "satellites2-40.mps",
-                                      "satellites2-60-fs.mps",
-                                      "savsched1.mps",
-                                      "sct2.mps",
-                                      "seymour.mps",
-                                      "seymour1.mps",
-                                      "sing326.mps",
-                                      "sing44.mps",
-                                      "sorrell3.mps",
-                                      "sp97ar.mps",
-                                      "sp98ar.mps",
-                                      "splice1k1.mps",
-                                      "square41.mps",
-                                      "square47.mps",
-                                      "supportcase10.mps",
-                                      "supportcase12.mps",
-                                      "supportcase18.mps",
-                                      "supportcase26.mps",
-                                      "supportcase33.mps",
-                                      "supportcase40.mps",
-                                      "supportcase42.mps",
-                                      "supportcase6.mps",
-                                      "supportcase7.mps",
-                                      "swath1.mps",
-                                      "swath3.mps",
-                                      "tbfp-network.mps",
-                                      "thor50dday.mps",
-                                      "timtab1.mps",
-                                      "tr12-30.mps",
-                                      "traininstance2.mps",
-                                      "traininstance6.mps",
-                                      "trento1.mps",
-                                      "uccase12.mps",
-                                      "uct-subprob.mps",
-                                      "unitcal_7.mps",
-                                      "var-smallemery-m6j6.mps"};
+std::vector<std::string> instances = {"supportcase26_presolved.mps",
+                                      "supportcase26_presolved_2.mps",
+                                      "supportcase26_presolved_3.mps",
+                                      "supportcase26_presolved_4.mps",
+                                      "supportcase26_presolved_5.mps",
+                                      "supportcase26_presolved_6.mps",
+                                      "supportcase26_presolved_7.mps",
+                                      "supportcase26_presolved_8.mps"};

benchmarks/linear_programming/cuopt/run_mip.cpp

@@ -210,6 +210,7 @@ int run_single_file(std::string file_path,
   settings.log_to_console                = log_to_console;
   settings.tolerances.relative_tolerance = 1e-12;
   settings.tolerances.absolute_tolerance = 1e-6;
+  settings.presolve                      = true;
   cuopt::linear_programming::benchmark_info_t benchmark_info;
   settings.benchmark_info_ptr = &benchmark_info;
   auto start_run_solver       = std::chrono::high_resolution_clock::now();

cpp/include/cuopt/linear_programming/mip/solver_settings.hpp

@@ -72,8 +72,8 @@ class mip_solver_settings_t {
 
   struct tolerances_t {
     f_t presolve_absolute_tolerance = 1.0e-6;
-    f_t absolute_tolerance          = 1.0e-4;
-    f_t relative_tolerance          = 1.0e-6;
+    f_t absolute_tolerance          = 1.0e-6;
+    f_t relative_tolerance          = 1.0e-12;
     f_t integrality_tolerance       = 1.0e-5;
     f_t absolute_mip_gap            = 1.0e-10;
     f_t relative_mip_gap            = 1.0e-4;

cpp/include/cuopt/logger.hpp

@@ -78,7 +78,7 @@ inline rapids_logger::logger& default_logger()
     logger_.set_pattern(default_pattern());
 #endif
     logger_.set_level(default_level());
-    logger_.flush_on(rapids_logger::level_enum::info);
+    logger_.flush_on(rapids_logger::level_enum::debug);
 
     return logger_;
   }();
@@ -100,7 +100,7 @@ inline void reset_default_logger()
   default_logger().set_pattern(default_pattern());
 #endif
   default_logger().set_level(default_level());
-  default_logger().flush_on(rapids_logger::level_enum::info);
+  default_logger().flush_on(rapids_logger::level_enum::debug);
 }
 
 }  // namespace cuopt

cpp/src/dual_simplex/bound_flipping_ratio_test.cpp

@@ -93,6 +93,11 @@ i_t bound_flipping_ratio_test_t<i_t, f_t>::single_pass(i_t start,
   }
   step_length       = min_val;
   nonbasic_entering = candidate;
+  // this should be temporary, find root causes where the candidate is not filled
+  if (nonbasic_entering == -1) {
+    // -1,-2 and -3 are reserved for other things
+    return -4;
+  }
   const i_t j = entering_index = nonbasic_list_[nonbasic_entering];
 
   constexpr bool verbose = false;
@@ -137,6 +142,7 @@ i_t bound_flipping_ratio_test_t<i_t, f_t>::compute_step_length(f_t& step_length,
 
   i_t k_idx = single_pass(
     0, num_breakpoints, indicies, ratios, slope, step_length, nonbasic_entering, entering_index);
+  if (k_idx == -4) { return -4; }
   bool continue_search = k_idx >= 0 && num_breakpoints > 1 && slope > 0.0;
   if (!continue_search) {
     if constexpr (0) {

cpp/src/dual_simplex/branch_and_bound.cpp

@@ -409,9 +409,12 @@ mip_status_t branch_and_bound_t<i_t, f_t>::solve(mip_solution_t<i_t, f_t>& solut
   assert(root_vstatus.size() == original_lp.num_cols);
   if (root_status == lp_status_t::INFEASIBLE) {
     settings.log.printf("MIP Infeasible\n");
-    if (settings.heuristic_preemption_callback != nullptr) {
-      settings.heuristic_preemption_callback();
-    }
+    // FIXME: rarely dual simplex detects infeasible whereas it is feasible.
+    // to add a small safety net, check if there is a primal solution already.
+    // Uncomment this if the issue with cost266-UUE is resolved
+    // if (settings.heuristic_preemption_callback != nullptr) {
+    //   settings.heuristic_preemption_callback();
+    // }
     return mip_status_t::INFEASIBLE;
   }
   if (root_status == lp_status_t::UNBOUNDED) {
@@ -434,8 +437,16 @@ mip_status_t branch_and_bound_t<i_t, f_t>::solve(mip_solution_t<i_t, f_t>& solut
   if (settings.set_simplex_solution_callback != nullptr) {
     std::vector<f_t> original_x;
     uncrush_primal_solution(original_problem, original_lp, root_relax_soln.x, original_x);
-    settings.set_simplex_solution_callback(original_x,
-                                           compute_user_objective(original_lp, root_objective));
+    std::vector<f_t> original_dual;
+    std::vector<f_t> original_z;
+    uncrush_dual_solution(original_problem,
+                          original_lp,
+                          root_relax_soln.y,
+                          root_relax_soln.z,
+                          original_dual,
+                          original_z);
+    settings.set_simplex_solution_callback(
+      original_x, original_dual, compute_user_objective(original_lp, root_objective));
   }
   mutex_lower.lock();
   f_t lower_bound = lower_bound_ = root_objective;

cpp/src/dual_simplex/phase2.cpp

@@ -1495,13 +1495,13 @@ void compute_delta_y(const basis_update_mpf_t<i_t, f_t>& ft,
 }
 
 template <typename i_t, typename f_t>
-void update_dual_variables(const sparse_vector_t<i_t, f_t>& delta_y_sparse,
-                           const std::vector<i_t>& delta_z_indices,
-                           const std::vector<f_t>& delta_z,
-                           f_t step_length,
-                           i_t leaving_index,
-                           std::vector<f_t>& y,
-                           std::vector<f_t>& z)
+i_t update_dual_variables(const sparse_vector_t<i_t, f_t>& delta_y_sparse,
+                          const std::vector<i_t>& delta_z_indices,
+                          const std::vector<f_t>& delta_z,
+                          f_t step_length,
+                          i_t leaving_index,
+                          std::vector<f_t>& y,
+                          std::vector<f_t>& z)
 {
   // Update dual variables
   // y <- y + steplength * delta_y
@@ -1517,6 +1517,7 @@ void update_dual_variables(const sparse_vector_t<i_t, f_t>& delta_y_sparse,
     z[j] += step_length * delta_z[j];
   }
   z[leaving_index] += step_length * delta_z[leaving_index];
+  return 0;
 }
 
 template <typename i_t, typename f_t>
@@ -2514,6 +2515,10 @@ dual::status_t dual_phase2(i_t phase,
                                                  delta_z_indices,
                                                  nonbasic_mark);
       entering_index = bfrt.compute_step_length(step_length, nonbasic_entering_index);
+      if (entering_index == -4) {
+        settings.log.printf("Numerical issues encountered in ratio test.\n");
+        return dual::status_t::NUMERICAL;
+      }
       timers.bfrt_time += timers.stop_timer();
     } else {
       entering_index = phase2::phase2_ratio_test(
@@ -2663,8 +2668,12 @@ dual::status_t dual_phase2(i_t phase,
     // Update dual variables
     // y <- y + steplength * delta_y
     // z <- z + steplength * delta_z
-    phase2::update_dual_variables(
+    i_t update_dual_variables_status = phase2::update_dual_variables(
       delta_y_sparse, delta_z_indices, delta_z, step_length, leaving_index, y, z);
+    if (update_dual_variables_status == -1) {
+      settings.log.printf("Numerical issues encountered in update_dual_variables.\n");
+      return dual::status_t::NUMERICAL;
+    }
     timers.vector_time += timers.stop_timer();
 
 #ifdef COMPUTE_DUAL_RESIDUAL

cpp/src/dual_simplex/presolve.cpp

@@ -1134,6 +1134,7 @@ void uncrush_dual_solution(const user_problem_t<i_t, f_t>& user_problem,
                            std::vector<f_t>& user_y,
                            std::vector<f_t>& user_z)
 {
+  user_y.resize(user_problem.num_rows);
   // Reduced costs are uncrushed just like the primal solution
   uncrush_primal_solution(user_problem, problem, z, user_z);
 

cpp/src/dual_simplex/simplex_solver_settings.hpp

@@ -109,7 +109,7 @@ struct simplex_solver_settings_t {
   i_t inside_mip;  // 0 if outside MIP, 1 if inside MIP at root node, 2 if inside MIP at leaf node
   std::function<void(std::vector<f_t>&, f_t)> solution_callback;
   std::function<void()> heuristic_preemption_callback;
-  std::function<void(std::vector<f_t>&, f_t)> set_simplex_solution_callback;
+  std::function<void(std::vector<f_t>&, std::vector<f_t>&, f_t)> set_simplex_solution_callback;
   mutable logger_t log;
   std::atomic<i_t>* concurrent_halt;  // if nullptr ignored, if !nullptr, 0 if solver should
                                       // continue, 1 if solver should halt

cpp/src/linear_programming/initial_scaling_strategy/initial_scaling.cu

@@ -39,17 +39,17 @@ pdlp_initial_scaling_strategy_t<i_t, f_t>::pdlp_initial_scaling_strategy_t(
   problem_t<i_t, f_t>& op_problem_scaled,
   i_t number_of_ruiz_iterations,
   f_t alpha,
-  pdhg_solver_t<i_t, f_t>& pdhg_solver,
   rmm::device_uvector<f_t>& A_T,
   rmm::device_uvector<i_t>& A_T_offsets,
   rmm::device_uvector<i_t>& A_T_indices,
+  pdhg_solver_t<i_t, f_t>* pdhg_solver_ptr,
   bool running_mip)
   : handle_ptr_(handle_ptr),
     stream_view_(handle_ptr_->get_stream()),
     primal_size_h_(op_problem_scaled.n_variables),
     dual_size_h_(op_problem_scaled.n_constraints),
     op_problem_scaled_(op_problem_scaled),
-    pdhg_solver_(pdhg_solver),
+    pdhg_solver_ptr_(pdhg_solver_ptr),
     A_T_(A_T),
     A_T_offsets_(A_T_offsets),
     A_T_indices_(A_T_indices),
@@ -398,7 +398,7 @@ void pdlp_initial_scaling_strategy_t<i_t, f_t>::scale_problem()
 
   op_problem_scaled_.is_scaled_ = true;
   if (!running_mip_) {
-    scale_solutions(pdhg_solver_.get_primal_solution(), pdhg_solver_.get_dual_solution());
+    scale_solutions(pdhg_solver_ptr_->get_primal_solution(), pdhg_solver_ptr_->get_dual_solution());
   }
 }