diff --git a/cpp/include/cuopt/linear_programming/constants.h b/cpp/include/cuopt/linear_programming/constants.h index ae02ee00c..024fae758 100644 --- a/cpp/include/cuopt/linear_programming/constants.h +++ b/cpp/include/cuopt/linear_programming/constants.h @@ -115,10 +115,12 @@ "mip_hyper_heuristic_related_vars_time_limit" /* @brief Diving heuristic toggles: -1 automatic, 0 disabled, 1 enabled */ -#define CUOPT_MIP_HYPER_DIVING_LINE_SEARCH "mip_hyper_diving_line_search" -#define CUOPT_MIP_HYPER_DIVING_PSEUDOCOST "mip_hyper_diving_pseudocost" -#define CUOPT_MIP_HYPER_DIVING_GUIDED "mip_hyper_diving_guided" -#define CUOPT_MIP_HYPER_DIVING_COEFFICIENT "mip_hyper_diving_coefficient" +#define CUOPT_MIP_HYPER_DIVING_LINE_SEARCH "mip_hyper_diving_line_search" +#define CUOPT_MIP_HYPER_DIVING_PSEUDOCOST "mip_hyper_diving_pseudocost" +#define CUOPT_MIP_HYPER_DIVING_GUIDED "mip_hyper_diving_guided" +#define CUOPT_MIP_HYPER_DIVING_COEFFICIENT "mip_hyper_diving_coefficient" +#define CUOPT_MIP_HYPER_DIVING_FARKAS "mip_hyper_diving_farkas" +#define CUOPT_MIP_HYPER_DIVING_VECTOR_LENGTH "mip_hyper_diving_vector_length" /* @brief Diving heuristic limits */ #define CUOPT_MIP_HYPER_DIVING_MIN_NODE_DEPTH "mip_hyper_diving_min_node_depth" #define CUOPT_MIP_HYPER_DIVING_NODE_LIMIT "mip_hyper_diving_node_limit" diff --git a/cpp/include/cuopt/linear_programming/mip/diving_hyper_params.hpp b/cpp/include/cuopt/linear_programming/mip/diving_hyper_params.hpp index 6ec005833..19463b915 100644 --- a/cpp/include/cuopt/linear_programming/mip/diving_hyper_params.hpp +++ b/cpp/include/cuopt/linear_programming/mip/diving_hyper_params.hpp @@ -20,10 +20,12 @@ namespace cuopt::linear_programming { template struct mip_diving_hyper_params_t { // -1 automatic, 0 disabled, 1 enabled - i_t line_search_diving = -1; - i_t pseudocost_diving = -1; - i_t guided_diving = -1; - i_t coefficient_diving = -1; + i_t line_search_diving = -1; + i_t pseudocost_diving = -1; + i_t guided_diving = -1; + i_t coefficient_diving = -1; + i_t farkas_diving = -1; + i_t vector_length_diving = -1; // The minimum depth to start diving from. i_t min_node_depth = 10; @@ -39,6 +41,11 @@ struct mip_diving_hyper_params_t { // The maximum backtracking allowed. i_t backtrack_limit = 5; + // For the Farkas diving to be effective, the coefficients in the objective function + // must have distinct values. The low tolerance here disables Farkas diving for + // set covering/partitioning, where all coefficients have the same value. + f_t farkas_obj_dynamism_tol = 1E-4; + // If a given diving heuristic found a new incumbent, show the corresponding // symbol in the first column of the log row. When false, every dive collapses // to 'D'. Otherwise, diff --git a/cpp/src/branch_and_bound/branch_and_bound.cpp b/cpp/src/branch_and_bound/branch_and_bound.cpp index a568aacd4..ac1a770b2 100644 --- a/cpp/src/branch_and_bound/branch_and_bound.cpp +++ b/cpp/src/branch_and_bound/branch_and_bound.cpp @@ -882,6 +882,12 @@ branch_variable_t branch_and_bound_t::variable_selection( mutex_upper_.unlock(); return guided_diving(pc_, fractional, solution, current_incumbent, log); + case FARKAS_DIVING: + return farkas_diving(worker->leaf_problem, fractional, solution, settings_.zero_tol, log); + + case VECTOR_LENGTH_DIVING: + return vector_length_diving(worker->leaf_problem, fractional, solution, log); + default: log.debug("Unknown variable selection method: %d\n", worker->search_strategy); return {-1, branch_direction_t::NONE}; @@ -1148,6 +1154,13 @@ struct deterministic_diving_policy_t log); } + case VECTOR_LENGTH_DIVING: + return vector_length_diving(this->worker.leaf_problem, fractional, x, log); + + case FARKAS_DIVING: + return farkas_diving( + this->worker.leaf_problem, fractional, x, this->bnb.settings_.zero_tol, log); + default: CUOPT_LOG_ERROR("Invalid diving method!"); return {-1, branch_direction_t::NONE}; } } @@ -1745,6 +1758,18 @@ void branch_and_bound_t::best_first_search_with(bfs_worker_t diving_settings.guided_diving = 0; } + if (diving_settings.farkas_diving != 0) { + f_t obj_dyn; + if (std::abs(original_lp_.min_abs_obj_coeff) < settings_.zero_tol) { + obj_dyn = std::abs(original_lp_.max_abs_obj_coeff) < settings_.zero_tol + ? 0 + : std::numeric_limits::infinity(); + } else { + obj_dyn = std::log10(original_lp_.max_abs_obj_coeff / original_lp_.min_abs_obj_coeff); + } + if (obj_dyn < diving_settings.farkas_obj_dynamism_tol) { diving_settings.farkas_diving = 0; } + } + worker->calculate_num_diving_workers( bfs_worker_pool_.size(), diving_worker_pool_.size(), diving_settings); diff --git a/cpp/src/branch_and_bound/constants.hpp b/cpp/src/branch_and_bound/constants.hpp index fed422ef8..d507999ed 100644 --- a/cpp/src/branch_and_bound/constants.hpp +++ b/cpp/src/branch_and_bound/constants.hpp @@ -9,23 +9,33 @@ namespace cuopt::linear_programming::dual_simplex { -constexpr int num_search_strategies = 5; +constexpr int num_search_strategies = 7; // Indicate the search and variable selection algorithms used by each thread // in B&B (See [1]). // // [1] T. Achterberg, “Constraint Integer Programming,” PhD, Technischen Universität Berlin, // Berlin, 2007. doi: 10.14279/depositonce-1634. +// [2] J. Witzig and A. Gleixner, “Conflict-Driven Heuristics for Mixed Integer Programming,” +// Feb. 07, 2019, _arXiv_: arXiv:1902.02615. doi: +// [10.48550/arXiv.1902.02615](https://doi.org/10.48550/arXiv.1902.02615). enum search_strategy_t : int { - BEST_FIRST = 0, // Best-First + Plunging. - PSEUDOCOST_DIVING = 1, // Pseudocost diving (9.2.5) - LINE_SEARCH_DIVING = 2, // Line search diving (9.2.4) - GUIDED_DIVING = 3, // Guided diving (9.2.3). - COEFFICIENT_DIVING = 4 // Coefficient diving (9.2.1) + BEST_FIRST = 0, // Best-First + Plunging. + PSEUDOCOST_DIVING = 1, // Pseudocost diving (9.2.5) + LINE_SEARCH_DIVING = 2, // Line search diving (9.2.4) + GUIDED_DIVING = 3, // Guided diving (9.2.3). + COEFFICIENT_DIVING = 4, // Coefficient diving (9.2.1) + FARKAS_DIVING = 5, // Farkas Diving (see [2]) + VECTOR_LENGTH_DIVING = 6 // Vector Length Diving (9.2.6) }; -constexpr search_strategy_t search_strategies[] = { - BEST_FIRST, PSEUDOCOST_DIVING, LINE_SEARCH_DIVING, GUIDED_DIVING, COEFFICIENT_DIVING}; +constexpr search_strategy_t search_strategies[] = {BEST_FIRST, + PSEUDOCOST_DIVING, + LINE_SEARCH_DIVING, + GUIDED_DIVING, + COEFFICIENT_DIVING, + FARKAS_DIVING, + VECTOR_LENGTH_DIVING}; enum class branch_direction_t { NONE = -1, DOWN = 0, UP = 1 }; diff --git a/cpp/src/branch_and_bound/diving_heuristics.cpp b/cpp/src/branch_and_bound/diving_heuristics.cpp index 46a681679..933f67f50 100644 --- a/cpp/src/branch_and_bound/diving_heuristics.cpp +++ b/cpp/src/branch_and_bound/diving_heuristics.cpp @@ -259,6 +259,112 @@ branch_variable_t coefficient_diving(const lp_problem_t& lp_probl return {branch_var, round_dir}; } +template +branch_variable_t farkas_diving(const lp_problem_t& lp, + const std::vector& fractional, + const std::vector& solution, + f_t zero_tol, + logger_t& log) +{ + if (fractional.size() == 0) return {-1, branch_direction_t::NONE}; + + i_t branch_var = -1; + f_t max_score = -1; + branch_direction_t round_dir = branch_direction_t::NONE; + + for (i_t j : fractional) { + f_t c = lp.objective[j]; + f_t f_down = solution[j] - std::floor(solution[j]); + f_t f_up = std::ceil(solution[j]) - solution[j]; + f_t score = 0; + branch_direction_t dir = branch_direction_t::NONE; + + if (c > zero_tol) { + dir = branch_direction_t::DOWN; + } else if (c < -zero_tol) { + dir = branch_direction_t::UP; + } else if (f_down < 0.5) { + dir = branch_direction_t::DOWN; + } else { + dir = branch_direction_t::UP; + } + + if (dir == branch_direction_t::UP) { + score = std::isfinite(lp.upper[j]) + ? (lp.upper[j] - std::floor(solution[j])) * f_up * std::abs(c) + : std::numeric_limits::infinity(); + } else { + score = std::isfinite(lp.lower[j]) + ? (std::ceil(solution[j]) - lp.lower[j]) * f_down * std::abs(c) + : std::numeric_limits::infinity(); + } + + if (score > max_score) { + max_score = score; + branch_var = j; + round_dir = dir; + } + } + + assert(round_dir != branch_direction_t::NONE); + assert(branch_var >= 0); + + log.debug("Farkas diving: selected var %d with val = %e, round dir = %d\n", + branch_var, + solution[branch_var], + round_dir); + + return {branch_var, round_dir}; +} + +template +branch_variable_t vector_length_diving(const lp_problem_t& lp, + const std::vector& fractional, + const std::vector& solution, + logger_t& log) +{ + if (fractional.size() == 0) return {-1, branch_direction_t::NONE}; + + constexpr f_t eps = 1E-6; + i_t branch_var = -1; + f_t min_score = std::numeric_limits::infinity(); + branch_direction_t round_dir = branch_direction_t::NONE; + + for (i_t j : fractional) { + f_t c = lp.objective[j]; + f_t f_down = solution[j] - std::floor(solution[j]); + f_t f_up = std::ceil(solution[j]) - solution[j]; + i_t column_length = lp.A.col_start[j + 1] - lp.A.col_start[j]; + f_t score = 0; + branch_direction_t dir = branch_direction_t::NONE; + + if (c < 0) { + dir = branch_direction_t::DOWN; + score = (f_down * std::abs(c) + eps) / (column_length + 1); + } else { + dir = branch_direction_t::UP; + score = (f_up * std::abs(c) + eps) / (column_length + 1); + } + + if (score < min_score) { + branch_var = j; + round_dir = dir; + min_score = score; + } + } + + assert(round_dir != branch_direction_t::NONE); + assert(branch_var >= 0); + + log.debug("Vector length diving: selected var %d with val = %e, round dir = %d and score = %e\n", + branch_var, + solution[branch_var], + round_dir, + min_score); + + return {branch_var, round_dir}; +} + #ifdef DUAL_SIMPLEX_INSTANTIATE_DOUBLE template branch_variable_t line_search_diving(const std::vector& fractional, const std::vector& solution, @@ -287,6 +393,18 @@ template branch_variable_t coefficient_diving(const lp_problem_t& up_locks, const std::vector& down_locks, logger_t& log); + +template branch_variable_t farkas_diving(const lp_problem_t& lp_problem, + const std::vector& fractional, + const std::vector& solution, + double zero_tol, + logger_t& log); + +template branch_variable_t vector_length_diving(const lp_problem_t& lp_problem, + const std::vector& fractional, + const std::vector& solution, + logger_t& log); + #endif } // namespace cuopt::linear_programming::dual_simplex diff --git a/cpp/src/branch_and_bound/diving_heuristics.hpp b/cpp/src/branch_and_bound/diving_heuristics.hpp index 79a23f08d..b2f1824cf 100644 --- a/cpp/src/branch_and_bound/diving_heuristics.hpp +++ b/cpp/src/branch_and_bound/diving_heuristics.hpp @@ -29,6 +29,8 @@ inline char feasible_solution_symbol(search_strategy_t strategy, bool log_diving case LINE_SEARCH_DIVING: return 'L'; case PSEUDOCOST_DIVING: return 'P'; case GUIDED_DIVING: return 'G'; + case FARKAS_DIVING: return 'F'; + case VECTOR_LENGTH_DIVING: return 'V'; default: return 'U'; } } @@ -43,6 +45,8 @@ bool is_search_strategy_enabled(search_strategy_t strategy, case LINE_SEARCH_DIVING: return settings.line_search_diving != 0; case GUIDED_DIVING: return settings.guided_diving != 0; case COEFFICIENT_DIVING: return settings.coefficient_diving != 0; + case FARKAS_DIVING: return settings.farkas_diving != 0; + case VECTOR_LENGTH_DIVING: return settings.vector_length_diving != 0; } return false; @@ -83,4 +87,17 @@ branch_variable_t coefficient_diving(const lp_problem_t& lp_probl const std::vector& down_locks, logger_t& log); +template +branch_variable_t farkas_diving(const lp_problem_t& lp, + const std::vector& fractional, + const std::vector& solution, + f_t zero_tol, + logger_t& log); + +template +branch_variable_t vector_length_diving(const lp_problem_t& lp, + const std::vector& fractional, + const std::vector& solution, + logger_t& log); + } // namespace cuopt::linear_programming::dual_simplex diff --git a/cpp/src/dual_simplex/presolve.hpp b/cpp/src/dual_simplex/presolve.hpp index c6587c9d7..22e578047 100644 --- a/cpp/src/dual_simplex/presolve.hpp +++ b/cpp/src/dual_simplex/presolve.hpp @@ -53,6 +53,11 @@ struct lp_problem_t { i_t cone_var_start{0}; std::vector second_order_cone_dims; + // Maximum and minimum value of the coefficients in the objective function. This is used + // for determine the "objective dynamism" in Farkas diving. + f_t max_abs_obj_coeff = 0; + f_t min_abs_obj_coeff = 0; + void write_mps(const std::string& path) const { std::ofstream mps_file(path); diff --git a/cpp/src/math_optimization/solver_settings.cu b/cpp/src/math_optimization/solver_settings.cu index df51dd688..4a000c6e0 100644 --- a/cpp/src/math_optimization/solver_settings.cu +++ b/cpp/src/math_optimization/solver_settings.cu @@ -168,9 +168,11 @@ solver_settings_t::solver_settings_t() : pdlp_settings(), mip_settings {CUOPT_MIP_HYPER_DIVING_PSEUDOCOST, &mip_settings.diving_params.pseudocost_diving, -1, 1, -1, "pseudocost diving toggle: -1 automatic, 0 disabled, 1 enabled"}, {CUOPT_MIP_HYPER_DIVING_GUIDED, &mip_settings.diving_params.guided_diving, -1, 1, -1, "guided diving toggle: -1 automatic, 0 disabled, 1 enabled"}, {CUOPT_MIP_HYPER_DIVING_COEFFICIENT, &mip_settings.diving_params.coefficient_diving, -1, 1, -1, "coefficient diving toggle: -1 automatic, 0 disabled, 1 enabled"}, + {CUOPT_MIP_HYPER_DIVING_FARKAS, &mip_settings.diving_params.farkas_diving, -1, 1, -1, "Farkas diving toggle: -1 automatic, 0 disabled, 1 enabled"}, + {CUOPT_MIP_HYPER_DIVING_VECTOR_LENGTH, &mip_settings.diving_params.vector_length_diving, -1, 1, -1, "vector-length diving toggle: -1 automatic, 0 disabled, 1 enabled"}, {CUOPT_MIP_HYPER_DIVING_MIN_NODE_DEPTH, &mip_settings.diving_params.min_node_depth, 0, std::numeric_limits::max(), 10, "minimum depth at which to start diving"}, {CUOPT_MIP_HYPER_DIVING_NODE_LIMIT, &mip_settings.diving_params.node_limit, 0, std::numeric_limits::max(), 500, "maximum nodes explored per dive"}, - {CUOPT_MIP_HYPER_DIVING_BACKTRACK_LIMIT, &mip_settings.diving_params.backtrack_limit, 0, std::numeric_limits::max(), 5, "maximum backtracking allowed per dive"}, + {CUOPT_MIP_HYPER_DIVING_BACKTRACK_LIMIT, &mip_settings.diving_params.backtrack_limit, 0, std::numeric_limits::max(), 5, "maximum backtracking allowed per dive"}, }; // Bool parameters