diff --git a/src/branch_and_cut.c b/src/branch_and_cut.c index 2e2ab3f..ef4baad 100644 --- a/src/branch_and_cut.c +++ b/src/branch_and_cut.c @@ -97,9 +97,7 @@ static int BNC_solve_node(struct BNC *bnc, int depth) rval = LP_get_obj_val(lp, &objval); abort_if(rval, "LP_get_obj_val failed\n"); - log_debug(" obj value = %.2f\n", objval); - - if (objval > *best_val + LP_EPSILON) + if (ceil(objval) > *best_val + LP_EPSILON) { log_debug("Branch pruned by bound (%.2lf > %.2lf).\n", objval, *best_val); @@ -120,20 +118,25 @@ static int BNC_solve_node(struct BNC *bnc, int depth) log_debug("Adding problem cutting planes...\n"); rval = bnc->problem_add_cutting_planes(lp, bnc->problem_data); abort_if(rval, "problem_add_cutting_planes failed"); - } - rval = LP_optimize(lp, &is_infeasible); - abort_if(rval, "LP_optimize failed\n"); + rval = LP_get_obj_val(lp, &objval); + abort_if(rval, "LP_get_obj_val failed"); - rval = LP_get_obj_val(lp, &objval); - abort_if(rval, "LP_get_obj_val failed"); + if (ceil(objval) > *best_val + LP_EPSILON) + { + log_debug("Branch pruned by bound (%.2lf > %.2lf).\n", objval, + *best_val); + rval = 0; + goto CLEANUP; + } - rval = LP_get_x(lp, x); - abort_if(rval, "LP_get_x failed"); + rval = LP_get_x(lp, x); + abort_if(rval, "LP_get_x failed"); + } if (BNC_is_integral(x, num_cols)) { - log_debug(" solution is integral\n"); + log_debug("Solution is integral\n"); if (objval + LP_EPSILON < *best_val) { @@ -148,9 +151,10 @@ static int BNC_solve_node(struct BNC *bnc, int depth) if (bnc->problem_solution_found) bnc->problem_solution_found(bnc->problem_data, bnc->best_x); } - } else + } + else { - log_debug(" solution is fractional\n"); + log_debug("Solution is fractional\n"); rval = BNC_branch_node(bnc, x, depth); abort_if(rval, "BNC_branch_node failed"); } diff --git a/src/gtsp-subtour.c b/src/gtsp-subtour.c index d9f9d28..fa9d1ea 100644 --- a/src/gtsp-subtour.c +++ b/src/gtsp-subtour.c @@ -172,9 +172,7 @@ int static add_subtour_cut( } int find_exact_subtour_cuts( - struct LP *lp, - struct GTSP *data, - double min_cut_violation) + struct LP *lp, struct GTSP *data, double min_cut_violation) { int rval = 0; @@ -215,8 +213,12 @@ int find_exact_subtour_cuts( abort_if(rval, "find_exact_subtour_cuts_cluster_to_cluster failed"); added_cuts_count = lp->cut_pool_size - original_cut_pool_size; - log_debug("Added %d cluster-to-cluster subtour cuts\n", added_cuts_count); - if (added_cuts_count > 0) goto CLEANUP; + if (added_cuts_count > 0) + { + log_debug("Added %d cluster-to-cluster subtour cuts\n", + added_cuts_count); + goto CLEANUP; + } // Constraints (2.2) original_cut_pool_size = lp->cut_pool_size; @@ -225,8 +227,11 @@ int find_exact_subtour_cuts( abort_if(rval, "find_exact_subtour_cuts_node_to_cluster failed"); added_cuts_count = lp->cut_pool_size - original_cut_pool_size; - log_debug("Added %d node-to-cluster subtour cuts\n", added_cuts_count); - if (added_cuts_count > 0) goto CLEANUP; + if (added_cuts_count > 0) + { + log_debug("Added %d node-to-cluster subtour cuts\n", added_cuts_count); + goto CLEANUP; + } // Constraints (2.3) original_cut_pool_size = lp->cut_pool_size; @@ -235,7 +240,11 @@ int find_exact_subtour_cuts( abort_if(rval, "find_exact_subtour_cuts_node_to_node failed"); added_cuts_count = lp->cut_pool_size - original_cut_pool_size; - log_debug("Added %d node-to-node subtour cuts\n", added_cuts_count); + if (added_cuts_count > 0) + { + log_debug("Added %d node-to-node subtour cuts\n", added_cuts_count); + goto CLEANUP; + } CLEANUP: graph_free(&digraph); diff --git a/src/gtsp.c b/src/gtsp.c index 27b53ba..c98ba56 100644 --- a/src/gtsp.c +++ b/src/gtsp.c @@ -213,34 +213,12 @@ int GTSP_add_cutting_planes(struct LP *lp, struct GTSP *data) } else { - log_debug("No more cuts found.\n"); + log_debug("No additional cuts found.\n"); break; } } - log_debug("Reoptimizing...\n"); - - double time_before_lp = get_current_time(); int is_infeasible; - rval = LP_optimize(lp, &is_infeasible); - abort_if(rval, "LP_optimize failed"); - double time_after_lp = get_current_time(); - - double obj_val; - rval = LP_get_obj_val(lp, &obj_val); - abort_if(rval, "LP_get_obj_val failed"); - - log_debug(" obj val = %.4lf\n", obj_val); - log_debug(" time = %.2lf\n", time_after_lp - time_before_lp); - - if (time_after_lp - time_before_lp > 10.0) - { - log_debug("LP is too slow. Removing slack constraints...\n"); - int start = data->graph->node_count + data->cluster_count; - rval = LP_remove_slacks(lp, start, LP_EPSILON); - abort_if(rval, "LP_remove_slacks failed"); - } - rval = LP_optimize(lp, &is_infeasible); abort_if(rval, "LP_optimize failed"); @@ -476,7 +454,8 @@ int GTSP_solution_found(struct GTSP *data, double *x) } double FLOW_CPU_TIME = 0; -double LP_CPU_TIME = 0; +double LP_SOLVE_TIME = 0; +double LP_CUT_POOL_TIME = 0; int LP_OPTIMIZE_COUNT = 0; int GTSP_main(int argc, char **argv) @@ -556,7 +535,8 @@ int GTSP_main(int argc, char **argv) log_info("Max-flow calls: %d\n", FLOW_MAX_FLOW_COUNT); log_info("Max-flow computation time: %.2lf\n", FLOW_CPU_TIME); log_info("LP optimize calls: %d\n", LP_OPTIMIZE_COUNT); - log_info("LP solving time: %.2lf\n", LP_CPU_TIME); + log_info("LP solving time: %.2lf\n", LP_SOLVE_TIME); + log_info("LP cut pool management time: %.2lf\n", LP_CUT_POOL_TIME); CLEANUP: GTSP_free(&data); diff --git a/src/lp.c b/src/lp.c index ad7bd9f..f2cd0dd 100644 --- a/src/lp.c +++ b/src/lp.c @@ -2,12 +2,10 @@ #include #include #include +#include #include "lp.h" #include "util.h" -#define LP_EPSILON 0.000001 -#define MAX_CUT_POOL_SIZE 100000 - int LP_open(struct LP *lp) { int rval = 0; @@ -16,6 +14,7 @@ int LP_open(struct LP *lp) lp->cplex_env = CPXopenCPLEX(&rval); abort_if(rval, "CPXopenCPLEX failed"); + lp->cut_pool_size = 0; lp->cut_pool = (struct Row **) malloc( MAX_CUT_POOL_SIZE * sizeof(struct Row *)); abort_if(!lp->cut_pool, "could not allocate cut_pool"); @@ -121,7 +120,42 @@ int LP_change_bound(struct LP *lp, int col, char lower_or_upper, double bnd) } extern int LP_OPTIMIZE_COUNT; -extern double LP_CPU_TIME; +extern double LP_SOLVE_TIME; +extern double LP_CUT_POOL_TIME; + +int LP_update_cut_ages(struct LP *lp) +{ + int rval = 0; + + double *slacks = 0; + + int numrows = CPXgetnumrows(lp->cplex_env, lp->cplex_lp); + + slacks = (double *) malloc(numrows * sizeof(double)); + abort_if(!slacks, "could not allocate slacks"); + + rval = CPXgetslack(lp->cplex_env, lp->cplex_lp, slacks, 0, numrows - 1); + abort_if(rval, "CPXgetslack failed"); + + int count = 0; + for (int i = 0; i < lp->cut_pool_size; i++) + { + struct Row *cut = lp->cut_pool[i]; + if (cut->cplex_row_index < 0) continue; + assert(cut->cplex_row_index < numrows); + + if (slacks[cut->cplex_row_index] < -LP_EPSILON) + cut->age++; + else + cut->age = 0; + + if (cut->age > 5) count++; + } + + CLEANUP: + if (slacks) free(slacks); + return rval; +} int LP_optimize(struct LP *lp, int *infeasible) { @@ -131,10 +165,17 @@ int LP_optimize(struct LP *lp, int *infeasible) *infeasible = 0; - double current = get_current_time(); + int numrows = CPXgetnumrows(lp->cplex_env, lp->cplex_lp); + int numcols = CPXgetnumcols(lp->cplex_env, lp->cplex_lp); + + log_debug("Optimizing LP (%d rows %d cols)...\n", numrows, numcols); + + double time_before = get_current_time(); rval = CPXdualopt(lp->cplex_env, lp->cplex_lp); abort_if(rval, "CPXdualopt failed"); - LP_CPU_TIME += get_current_time() - current; + + double time_after = get_current_time(); + LP_SOLVE_TIME += time_after - time_before; solstat = CPXgetstat(lp->cplex_env, lp->cplex_lp); if (solstat == CPX_STAT_INFEASIBLE) @@ -148,38 +189,82 @@ int LP_optimize(struct LP *lp, int *infeasible) solstat != CPX_STAT_OPTIMAL_INFEAS, "Invalid solution status"); } + double objval; + rval = LP_get_obj_val(lp, &objval); + abort_if(rval, "LP_get_obj_val failed"); + + log_debug(" obj val = %.4lf\n", objval); + log_debug(" time = %.4lf\n", time_after - time_before); + + time_before = get_current_time(); + rval = LP_update_cut_ages(lp); + abort_if(rval, "LP_update_cut_ages failed"); + + rval = LP_remove_old_cuts(lp); + abort_if(rval, "LP_remove_old_cuts failed"); + + time_after = get_current_time(); + LP_CUT_POOL_TIME += time_after - time_before; + CLEANUP: return rval; } -int LP_remove_slacks(struct LP *lp, int first_row, double max_slack) +int LP_remove_old_cuts(struct LP *lp) { int rval = 0; - double *slacks = 0; int *should_remove = 0; int numrows = CPXgetnumrows(lp->cplex_env, lp->cplex_lp); - if (numrows < 5000) return 0; + log_verbose(" numrows=%d\n", numrows); should_remove = (int *) malloc((numrows + 1) * sizeof(int)); abort_if(!should_remove, "could not allocate should_remove"); - slacks = (double *) malloc(numrows * sizeof(double)); - abort_if(!slacks, "could not allocate slacks"); + for (int i = 0; i < numrows; i++) + should_remove[i] = 0; + should_remove[numrows] = 0; - rval = CPXgetslack(lp->cplex_env, lp->cplex_lp, slacks, 0, numrows - 1); - abort_if(rval, "CPXgetslack failed"); + log_verbose("Old cplex row index:\n"); + for (int i = 0; i < lp->cut_pool_size; i++) + log_verbose(" %d\n", lp->cut_pool[i]->cplex_row_index); + + log_verbose("Should remove:\n"); + for (int i = 0; i < lp->cut_pool_size; i++) + { + struct Row *cut = lp->cut_pool[i]; + if (cut->age <= MAX_CUT_AGE) continue; + if (cut->cplex_row_index < 0) continue; + should_remove[cut->cplex_row_index] = 1; + log_verbose(" %d\n", cut->cplex_row_index); + } + + int count = 0; for (int i = 0; i < numrows; i++) - should_remove[i] = (slacks[i] < -max_slack); - should_remove[numrows] = 0; + { + if (!should_remove[i]) continue; + + for (int j = 0; j < lp->cut_pool_size; j++) + { + struct Row *cut = lp->cut_pool[j]; + + if (cut->cplex_row_index == i - count) cut->cplex_row_index = 0; + else if (cut->cplex_row_index > i - count) cut->cplex_row_index--; + } + + count++; + } + + log_verbose("New cplex row index:\n"); + for (int i = 0; i < lp->cut_pool_size; i++) + log_verbose(" %d\n", lp->cut_pool[i]->cplex_row_index); - log_debug("Deleting constraints...\n"); int start = 0; int end = -1; - int count = 0; - for (int i = first_row; i < numrows; i++) + count = 0; + for (int i = 0; i < numrows + 1; i++) { if (should_remove[i]) { @@ -202,11 +287,24 @@ int LP_remove_slacks(struct LP *lp, int first_row, double max_slack) } } - log_info("Removed %d of %d constraints\n", count, numrows); + numrows = CPXgetnumrows(lp->cplex_env, lp->cplex_lp); + log_verbose(" new numrows=%d\n", numrows); + + for (int i = 0; i < lp->cut_pool_size; i++) + { + struct Row *cut = lp->cut_pool[i]; + assert(cut->cplex_row_index < numrows); + } + + if (count > 0) + { + log_info("Found and removed %d old cuts\n", count); + rval = CPXdualopt(lp->cplex_env, lp->cplex_lp); + abort_if(rval, "CPXoptimize failed"); + } CLEANUP: if (should_remove) free(should_remove); - if (slacks) free(slacks); return rval; } @@ -280,11 +378,13 @@ int LP_add_cut(struct LP *lp, struct Row *cut) { int rval = 0; + double time_before = get_current_time(); rval = LP_update_hash(cut); abort_if(rval, "LP_update_hash failed"); for (int i = 0; i < lp->cut_pool_size; i++) { + if (lp->cut_pool[i]->cplex_row_index < 0) continue; if (lp->cut_pool[i]->hash != cut->hash) continue; if (!compare_cuts(lp->cut_pool[i], cut)) { @@ -295,6 +395,7 @@ int LP_add_cut(struct LP *lp, struct Row *cut) } } + abort_if(lp->cut_pool_size > MAX_CUT_POOL_SIZE, "Cut pool is too large"); lp->cut_pool[lp->cut_pool_size++] = cut; int rmatbeg = 0; @@ -302,6 +403,12 @@ int LP_add_cut(struct LP *lp, struct Row *cut) cut->rmatind, cut->rmatval); abort_if(rval, "LP_add_rows failed"); + cut->cplex_row_index = CPXgetnumrows(lp->cplex_env, lp->cplex_lp) - 1; + cut->age = 0; + + double time_after = get_current_time(); + LP_CUT_POOL_TIME += time_after - time_before; + CLEANUP: return rval; } diff --git a/src/lp.h b/src/lp.h index 64e35f5..de5f0d0 100644 --- a/src/lp.h +++ b/src/lp.h @@ -4,6 +4,8 @@ #include #define LP_EPSILON 0.000001 +#define MAX_CUT_POOL_SIZE 1000000 +#define MAX_CUT_AGE 10 struct LP { @@ -17,7 +19,8 @@ struct LP struct Row { unsigned long hash; -// int cplex_row_index; + int cplex_row_index; + int age; int nz; char sense; @@ -69,7 +72,7 @@ int LP_get_x(struct LP *lp, double *x); int LP_get_num_cols(struct LP *lp); -int LP_remove_slacks(struct LP *lp, int start, double max_slack); +int LP_remove_old_cuts(struct LP *lp); int LP_update_hash(struct Row *cut);