Skip to content

Commit

Permalink
Merge pull request #1989 from ERGO-Code/fix-1865
Browse files Browse the repository at this point in the history
Phase 2 primal simplex being forced inadvisedly
  • Loading branch information
jajhall authored Oct 27, 2024
2 parents 4d0949a + 6d3fbe2 commit a755161
Show file tree
Hide file tree
Showing 3 changed files with 70 additions and 27 deletions.
80 changes: 57 additions & 23 deletions src/simplex/HApp.h
Original file line number Diff line number Diff line change
Expand Up @@ -254,12 +254,27 @@ inline HighsStatus solveLpSimplex(HighsLpSolverObject& solver_object) {
scaled_model_status == HighsModelStatus::kOptimal &&
(num_unscaled_primal_infeasibilities ||
num_unscaled_dual_infeasibilities);
if (scaled_optimality_but_unscaled_infeasibilities)
// Determine whether the unscaled solution has primal
// infeasibilities after the scaled LP has been solved to the
// objective target
const bool scaled_objective_target_but_unscaled_primal_infeasibilities =
scaled_model_status == HighsModelStatus::kObjectiveTarget &&
highs_info.num_primal_infeasibilities > 0;
// Determine whether the unscaled solution has dual
// infeasibilities after the scaled LP has been solved to the
// objective bound
const bool scaled_objective_bound_but_unscaled_dual_infeasibilities =
scaled_model_status == HighsModelStatus::kObjectiveBound &&
highs_info.num_dual_infeasibilities > 0;
if (scaled_optimality_but_unscaled_infeasibilities ||
scaled_objective_target_but_unscaled_primal_infeasibilities ||
scaled_objective_bound_but_unscaled_dual_infeasibilities)
highsLogDev(options.log_options, HighsLogType::kInfo,
"Have num/max/sum primal (%" HIGHSINT_FORMAT
"/%g/%g) and dual (%" HIGHSINT_FORMAT
"After unscaling with status %s, have num/max/sum primal "
"(%" HIGHSINT_FORMAT "/%g/%g) and dual (%" HIGHSINT_FORMAT
"/%g/%g) "
"unscaled infeasibilities\n",
utilModelStatusToString(scaled_model_status).c_str(),
highs_info.num_primal_infeasibilities,
highs_info.max_primal_infeasibility,
highs_info.sum_primal_infeasibilities,
Expand All @@ -274,8 +289,8 @@ inline HighsStatus solveLpSimplex(HighsLpSolverObject& solver_object) {
scaled_model_status == HighsModelStatus::kInfeasible ||
scaled_model_status == HighsModelStatus::kUnboundedOrInfeasible ||
scaled_model_status == HighsModelStatus::kUnbounded ||
scaled_model_status == HighsModelStatus::kObjectiveBound ||
scaled_model_status == HighsModelStatus::kObjectiveTarget ||
scaled_objective_bound_but_unscaled_dual_infeasibilities ||
scaled_objective_target_but_unscaled_primal_infeasibilities ||
scaled_model_status == HighsModelStatus::kUnknown);
// Handle the case when refinement will not take place
if (!refine_solution) {
Expand Down Expand Up @@ -322,21 +337,31 @@ inline HighsStatus solveLpSimplex(HighsLpSolverObject& solver_object) {
options.dual_simplex_cost_perturbation_multiplier;
HighsInt simplex_dual_edge_weight_strategy =
ekk_info.dual_edge_weight_strategy;
// #1865 exposed that this should not be
// HighsModelStatus::kObjectiveBound, but
// HighsModelStatus::kObjectiveTarget, since if the latter is
// the model status for the scaled LP, any primal
// infeasibilities should be small, but must be cleaned up
// before (hopefully) a few phase 2 primal simplex iterations
// are required to attain the target for the unscaled LP
//
// In #1865, phase 2 primal simplex was forced with large primal
// infeasibilities
if (num_unscaled_primal_infeasibilities == 0 ||
scaled_model_status == HighsModelStatus::kObjectiveBound) {
// Only dual infeasibilities, or primal infeasibilities do not
// matter due to solution status, so use primal simplex phase
// 2
scaled_model_status == HighsModelStatus::kObjectiveTarget) {
// Only dual infeasibilities, or objective target reached (in
// primal phase 2) - so primal infeasibilities should be small
// - so use primal simplex phase 2
options.simplex_strategy = kSimplexStrategyPrimal;
if (scaled_model_status == HighsModelStatus::kObjectiveBound) {
highsLogDev(
options.log_options, HighsLogType::kInfo,
"solveLpSimplex: Calling primal simplex after "
"scaled_model_status == HighsModelStatus::kObjectiveBound: solve "
"= %d; tick = %d; iter = %d\n",
(int)ekk_instance.debug_solve_call_num_,
(int)ekk_instance.debug_initial_build_synthetic_tick_,
(int)ekk_instance.iteration_count_);
if (scaled_model_status == HighsModelStatus::kObjectiveTarget) {
highsLogDev(options.log_options, HighsLogType::kInfo,
"solveLpSimplex: Calling primal simplex after "
"scaled_model_status == "
"HighsModelStatus::kObjectiveTarget: solve "
"= %d; tick = %d; iter = %d\n",
(int)ekk_instance.debug_solve_call_num_,
(int)ekk_instance.debug_initial_build_synthetic_tick_,
(int)ekk_instance.iteration_count_);
}
} else {
// Using dual simplex, so force Devex if starting from an advanced
Expand All @@ -349,11 +374,21 @@ inline HighsStatus solveLpSimplex(HighsLpSolverObject& solver_object) {
//
// Solve the unscaled LP with scaled NLA
//
// Force the simplex solver to start in phase 2 unless solving
// the LP directly as unscaled
// Force the simplex solver to start in phase 1 if solving the
// LP directly as unscaled, or using primal simplex to clean up
// small dual infeasibilities after the scaled LP yielded model
// status HighsModelStatus::kObjectiveTarget. Otherwise force
// the simplex solver to start in phase 2
//
const bool force_phase2 = options.simplex_unscaled_solution_strategy !=
kSimplexUnscaledSolutionStrategyDirect;
const bool force_phase1 =
(options.simplex_unscaled_solution_strategy ==
kSimplexUnscaledSolutionStrategyDirect) ||
(scaled_model_status == HighsModelStatus::kObjectiveTarget);
const bool force_phase2 =
(options.simplex_unscaled_solution_strategy !=
kSimplexUnscaledSolutionStrategyDirect) &&
(scaled_model_status != HighsModelStatus::kObjectiveTarget);
assert(force_phase2 == !force_phase1);
return_status = ekk_instance.solve(force_phase2);
solved_unscaled_lp = true;
if (scaled_model_status != HighsModelStatus::kObjectiveBound &&
Expand All @@ -362,7 +397,6 @@ inline HighsStatus solveLpSimplex(HighsLpSolverObject& solver_object) {
// for the first time in which case we again call solve with primal
// simplex if not dual feasible
const bool objective_bound_refinement =
ekk_instance.model_status_ == HighsModelStatus::kObjectiveBound &&
ekk_info.num_dual_infeasibilities > 0;
if (objective_bound_refinement) {
options.simplex_strategy = kSimplexStrategyPrimal;
Expand Down
8 changes: 4 additions & 4 deletions src/simplex/HEkkDual.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2807,10 +2807,10 @@ bool HEkkDual::reachedExactObjectiveBound() {
exact_dual_objective_value - objective_bound;
std::string action;
if (exact_dual_objective_value > objective_bound) {
highsLogDev(ekk_instance_.options_->log_options, HighsLogType::kDetailed,
"HEkkDual::solvePhase2: %12g = Objective > ObjectiveUB\n",
ekk_instance_.info_.updated_dual_objective_value,
objective_bound);
highsLogDev(
ekk_instance_.options_->log_options, HighsLogType::kDetailed,
"HEkkDual::solvePhase2: %12g = Objective > ObjectiveUB = %12g\n",
ekk_instance_.info_.updated_dual_objective_value, objective_bound);
action = "Have DualUB bailout";
if (ekk_instance_.info_.costs_perturbed ||
ekk_instance_.info_.costs_shifted) {
Expand Down
9 changes: 9 additions & 0 deletions src/simplex/HEkkPrimal.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2812,6 +2812,15 @@ void HEkkPrimal::shiftBound(const bool lower, const HighsInt iVar,
shift = infeasibility + feasibility;
bound -= shift;
new_infeasibility = bound - value;
if (new_infeasibility >= 0) {
printf(
"HEkkPrimal::shiftBound LB = %g; random_value = %g; value = %g; "
"feasibility = %g; infeasibility = %g; shift = %g; bound = %g; "
"new_infeasibility = %g; \n",
old_bound, random_value, value, feasibility, infeasibility, shift,
bound, new_infeasibility);
fflush(stdout);
}
assert(new_infeasibility < 0);
} else {
// Bound to shift is upper
Expand Down

0 comments on commit a755161

Please sign in to comment.