isoron
/
gtsp
Archived
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Simplify LP; remove useless cutting planes

Browse Source
master
Alinson S. Xavier 11 years ago
parent 1711695e19
commit b802256426
4 changed files with 44 additions and 294 deletions
Show Stats Download Patch File Download Diff File

27
src/gtsp-comb.c
Unescape View File

@ -41,7 +41,7 @@ int add_comb_cut(
if (components[clusters[e->from->index]] != current_component) continue;
if (components[clusters[e->to->index]] != current_component) continue;
rmatind[nz] = node_count + e->index;
rmatind[nz] = e->index;
rmatval[nz] = -1.0;
nz++;
@ -62,32 +62,11 @@ int add_comb_cut(
log_verbose(" tooth (%d %d)\n", e->from->index, e->to->index);
rmatind[nz] = node_count + e->index;
rmatind[nz] = e->index;
rmatval[nz] = -1.0;
nz++;
}
// // Lifting of the nodes
// for (int i = 0; i < node_count; i++)
// {
// double val;
// struct Node *n = &graph->nodes[i];
// int c = node_to_cluster[n->index];
//
// if (components[c] == current_component)
// val = (teeth[c] < 0 ? 1.0 : 0.0);
// else
// val = (teeth[c] < 0 ? 0.0 : 0.0);
//
// if (val == 0.0) continue;
//
// rmatind[nz] = n->index;
// rmatval[nz] = val;
// nz++;
//
// rhs = val;
// }
#if LOG_LEVEL >= LOG_LEVEL_DEBUG
log_verbose("Generated cut:\n");
if (OPTIMAL_X)
@ -357,7 +336,7 @@ static int shrink_clusters(
int to = clusters[e->to->index];
int shunk_e_index = edge_map[from * cluster_count + to];
shrunken_x[shunk_e_index] += x[graph->node_count + e->index];
shrunken_x[shunk_e_index] += x[e->index];
}
CLEANUP:

229
src/gtsp-subtour.c
Unescape View File

@ -28,7 +28,7 @@ int static build_flow_digraph(
int kc = 0;
for (int i = 0; i < graph->edge_count; i++)
{
if (x[node_count + i] < LP_EPSILON) continue;
if (x[i] < LP_EPSILON) continue;
struct Edge *e = &graph->edges[i];
int from = e->from->index;
@ -36,7 +36,7 @@ int static build_flow_digraph(
digraph_edges[ke++] = from;
digraph_edges[ke++] = to;
capacities[kc++] = x[node_count + i];
capacities[kc++] = x[i];
digraph_edges[ke++] = to;
digraph_edges[ke++] = from;
@ -44,7 +44,7 @@ int static build_flow_digraph(
digraph_edges[ke++] = to;
digraph_edges[ke++] = from;
capacities[kc++] = x[node_count + i];
capacities[kc++] = x[i];
digraph_edges[ke++] = from;
digraph_edges[ke++] = to;
@ -99,20 +99,16 @@ int static build_flow_digraph(
return rval;
}
int static add_subtour_cut(
int add_subtour_cut(
struct LP *lp,
struct Graph *graph,
struct Node *from,
struct Node *to,
struct Edge **cut_edges,
int cut_edges_count,
int type)
int cut_edges_count)
{
int rval = 0;
char sense = 'G';
double rhs = 2.0 - 2.0 * type;
int newnz = cut_edges_count + type;
double rhs = 2.0;
int newnz = cut_edges_count;
int *rmatind = 0;
double *rmatval = 0;
@ -125,22 +121,10 @@ int static add_subtour_cut(
for (int i = 0; i < cut_edges_count; i++)
{
rmatind[i] = cut_edges[i]->index + graph->node_count;
rmatind[i] = cut_edges[i]->index;
rmatval[i] = 1.0;
}
if (type >= 1)
{
rmatind[cut_edges_count] = from->index;
rmatval[cut_edges_count] = -2.0;
}
if (type >= 2)
{
rmatind[cut_edges_count + 1] = to->index;
rmatval[cut_edges_count + 1] = -2.0;
}
log_verbose("Generated cut:\n");
for (int i = 0; i < newnz; i++)
log_verbose(" %8.2f x%d\n", rmatval[i], rmatind[i]);
@ -208,7 +192,6 @@ int find_exact_subtour_cuts(
rval = build_flow_digraph(data, x, &digraph, capacities);
abort_if(rval, "build_flow_digraph failed");
// Constraints (2.1)
rval = find_exact_subtour_cuts_cluster_to_cluster(lp, data, &digraph,
capacities, min_cut_violation);
abort_if(rval, "find_exact_subtour_cuts_cluster_to_cluster failed");
@ -219,30 +202,6 @@ int find_exact_subtour_cuts(
if (added_cuts_count > 0)
goto CLEANUP;
// Constraints (2.2)
original_cut_pool_size = lp->cut_pool_size;
rval = find_exact_subtour_cuts_node_to_cluster(lp, data, x, &digraph,
capacities, min_cut_violation);
abort_if(rval, "find_exact_subtour_cuts_node_to_cluster failed");
added_cuts_count = lp->cut_pool_size - original_cut_pool_size;
log_debug(" %d node-to-cluster\n", added_cuts_count);
SUBTOUR_NODE_CLUSTER_COUNT += added_cuts_count;
if (added_cuts_count > 0)
goto CLEANUP;
// Constraints (2.3)
original_cut_pool_size = lp->cut_pool_size;
rval = find_exact_subtour_cuts_node_to_node(lp, data, x, &digraph,
capacities, min_cut_violation);
abort_if(rval, "find_exact_subtour_cuts_node_to_node failed");
added_cuts_count = lp->cut_pool_size - original_cut_pool_size;
log_debug(" %d node-to-node\n", added_cuts_count);
SUBTOUR_NODE_NODE_COUNT += added_cuts_count;
if (added_cuts_count > 0)
goto CLEANUP;
SUBTOUR_TIME += get_current_time() - initial_time;
CLEANUP:
@ -252,175 +211,6 @@ int find_exact_subtour_cuts(
return rval;
}
int find_exact_subtour_cuts_node_to_node(
struct LP *lp,
struct GTSP *data,
double *x,
struct Graph *digraph,
double *capacities,
double min_cut_violation)
{
int rval = 0;
struct Edge **cut_edges = 0;
double *flow = 0;
struct Graph *graph = data->graph;
int *clusters = data->node_to_cluster;
cut_edges = (struct Edge **) malloc(
graph->edge_count * sizeof(struct Edge *));
flow = (double *) malloc(digraph->edge_count * sizeof(double));
abort_if(!cut_edges, "could not allocate cut_edges");
abort_if(!flow, "could not allocate flow");
int max_x_index = 0;
double max_x = DBL_MIN;
for (int i = 0; i < graph->node_count; i++)
{
struct Node *n = &graph->nodes[i];
if (x[n->index] > max_x)
{
max_x = x[n->index];
max_x_index = i;
}
}
int i = max_x_index;
for (int j = 0; j < graph->node_count; j++)
{
if (i == j) continue;
if (clusters[i] == clusters[j]) continue;
if (x[i] + x[j] - 1 <= LP_EPSILON) continue;
struct Node *from = &digraph->nodes[i];
struct Node *to = &digraph->nodes[j];
int cut_edges_count;
double flow_value;
rval = flow_find_max_flow(digraph, capacities, from, to, flow,
&flow_value);
abort_if(rval, "flow_find_max_flow failed");
if (flow_value >= 2 * (x[i] + x[j] - 1) - min_cut_violation)
continue;
log_verbose("Marked nodes:\n");
for (int k = 0; k < graph->node_count; k++)
{
graph->nodes[k].mark = digraph->nodes[k].mark;
if (digraph->nodes[k].mark) log_verbose(" %d\n", k);
}
rval = get_cut_edges_from_marks(graph, &cut_edges_count, cut_edges);
abort_if(rval, "get_cut_edges_from_marks failed");
log_verbose("Cut edges:\n");
for (int k = 0; k < cut_edges_count; k++)
log_verbose(" %d %d (%d)\n", cut_edges[k]->from->index,
cut_edges[k]->to->index, cut_edges[k]->index);
rval = add_subtour_cut(lp, graph, from, to, cut_edges, cut_edges_count,
2);
abort_if(rval, "add_subtour_cut failed");
}
CLEANUP:
if (flow) free(flow);
if (cut_edges) free(cut_edges);
return rval;
}
int find_exact_subtour_cuts_node_to_cluster(
struct LP *lp,
struct GTSP *data,
double *x,
struct Graph *digraph,
double *capacities,
double min_cut_violation)
{
int rval = 0;
int cuts_count = 0;
struct Edge **cut_edges = 0;
double *flow = 0;
struct Graph *graph = data->graph;
int *clusters = data->node_to_cluster;
cut_edges = (struct Edge **) malloc(
graph->edge_count * sizeof(struct Edge *));
flow = (double *) malloc(digraph->edge_count * sizeof(double));
abort_if(!cut_edges, "could not allocate cut_edges");
abort_if(!flow, "could not allocate flow");
for (int i = 0; i < graph->node_count; i++)
{
for (int j = 0; j < data->cluster_count; j++)
{
if (clusters[i] == j) continue;
if (x[i] < LP_EPSILON) continue;
struct Node *from = &digraph->nodes[i];
struct Node *to = &digraph->nodes[graph->node_count + j];
log_verbose(
"Sending flow from node %d to cluster %d (must be >= %.4lf)\n",
i, j, 2 * x[i]);
activate_cluster_node(capacities, to);
double flow_value;
int cut_edges_count;
rval = flow_find_max_flow(digraph, capacities, from, to, flow,
&flow_value);
abort_if(rval, "flow_find_max_flow failed");
log_verbose(" flow value = %.4lf\n", flow_value);
deactivate_cluster_node(capacities, to);
if (flow_value + min_cut_violation >= 2 * x[i])
continue;
log_verbose("Marked nodes:\n");
for (int k = 0; k < graph->node_count; k++)
{
graph->nodes[k].mark = digraph->nodes[k].mark;
if (graph->nodes[k].mark) log_verbose(" %d\n", k);
}
rval = get_cut_edges_from_marks(graph, &cut_edges_count, cut_edges);
abort_if(rval, "get_cut_edges_from_marks failed");
log_verbose("Cut edges:\n");
for (int k = 0; k < cut_edges_count; k++)
{
struct Edge *e = cut_edges[k];
assert(e->from->mark != e->to->mark);
log_verbose(" %d (%d) %d (%d) [%d]\n", e->from->index,
e->from->mark, e->to->index, e->to->mark, e->index);
}
rval = add_subtour_cut(lp, graph, from, 0, cut_edges,
cut_edges_count, 1);
abort_if(rval, "add_subtour_cut failed");
cuts_count++;
}
}
CLEANUP:
if (cut_edges) free(cut_edges);
if (flow) free(flow);
return rval;
}
int find_exact_subtour_cuts_cluster_to_cluster(
struct LP *lp,
struct GTSP *data,
@ -489,8 +279,7 @@ int find_exact_subtour_cuts_cluster_to_cluster(
log_verbose(" %d %d (%d)\n", cut_edges[k]->from->index,
cut_edges[k]->to->index, cut_edges[k]->index);
rval = add_subtour_cut(lp, graph, 0, 0, cut_edges, cut_edges_count,
0);
rval = add_subtour_cut(lp, cut_edges, cut_edges_count);
abort_if(rval, "add_subtour_cut failed");
cuts_count++;

80
src/gtsp.c
Unescape View File

@ -4,6 +4,7 @@
#include <math.h>
#include <assert.h>
#include <sys/stat.h>
#include <float.h>
#include "gtsp.h"
#include "geometry.h"
#include "util.h"
@ -180,21 +181,14 @@ int GTSP_init_lp(struct LP *lp, struct GTSP *data)
{
int rval = 0;
int node_count = data->graph->node_count;
int edge_count = data->graph->edge_count;
int cluster_count = data->cluster_count;
int *clusters = data->node_to_cluster;
struct Edge *edges = data->graph->edges;
for (int i = 0; i < node_count; i++)
{
rval = LP_new_row(lp, 'E', 0.0);
abort_if(rval, "LP_new_row failed");
}
for (int i = 0; i < cluster_count; i++)
{
rval = LP_new_row(lp, 'E', 1.0);
rval = LP_new_row(lp, 'E', 2.0);
abort_if(rval, "LP_new_row failed");
}
@ -202,22 +196,14 @@ int GTSP_init_lp(struct LP *lp, struct GTSP *data)
double ub = 1.0;
int cmatbeg = 0;
for (int i = 0; i < node_count; i++)
{
double obj = 0.0;
double cmatval[] = {-2.0, 1.0};
int cmatind[] = {i, node_count + clusters[i]};
rval = LP_add_cols(lp, 1, 2, &obj, &cmatbeg, cmatind, cmatval, &lb,
&ub);
abort_if(rval, "LP_add_cols failed");
}
for (int i = 0; i < edge_count; i++)
{
struct Node *from = edges[i].from;
struct Node *to = edges[i].to;
double obj = (double) edges[i].weight;
double cmatval[] = {1.0, 1.0};
int cmatind[] = {edges[i].from->index, edges[i].to->index};
int cmatind[] = {clusters[from->index], clusters[to->index]};
rval = LP_add_cols(lp, 1, 2, &obj, &cmatbeg, cmatind, cmatval, &lb,
&ub);
@ -317,7 +303,7 @@ int GTSP_write_solution(struct GTSP *data, char *filename, double *x)
int positive_edge_count = 0;
for (int i = 0; i < edge_count; i++)
if (x[i + node_count] > LP_EPSILON)
if (x[i] > LP_EPSILON)
positive_edge_count++;
fprintf(file, "%d %d\n", node_count, edge_count);
@ -325,9 +311,9 @@ int GTSP_write_solution(struct GTSP *data, char *filename, double *x)
fprintf(file, "%d\n", positive_edge_count);
for (int i = 0; i < edge_count; i++)
if (x[i + node_count] > LP_EPSILON)
if (x[i] > LP_EPSILON)
fprintf(file, "%d %d %.4lf\n", edges[i].from->index,
edges[i].to->index, x[i + node_count]);
edges[i].to->index, x[i]);
CLEANUP:
if (file) fclose(file);
@ -354,12 +340,10 @@ int GTSP_read_solution(struct GTSP *gtsp, char *filename, double **p_x)
rval = fscanf(file, "%d %d", &node_count, &edge_count);
abort_if(rval != 2, "invalid input format (node and edge count)");
int num_cols = node_count + edge_count;
x = (double *) malloc(num_cols * sizeof(double));
x = (double *) malloc(edge_count * sizeof(double));
abort_if(!x, "could not allocate x");
for (int i = 0; i < node_count + edge_count; i++) x[i] = 0.0;
for (int i = 0; i < edge_count; i++) x[i] = 0.0;
rval = fscanf(file, "%d", &edge_count);
abort_if(rval != 1, "invalid input format (positive edge count)");
@ -378,14 +362,11 @@ int GTSP_read_solution(struct GTSP *gtsp, char *filename, double **p_x)
abort_if(rval != 3, "invalid input format (edge endpoints)");
edge = edge_map[from * node_count + to];
abort_if(edge > num_cols, "invalid edge");
x[from] += val / 2;
x[to] += val / 2;
x[edge] = val;
}
for (int i = 0; i < num_cols; i++)
for (int i = 0; i < edge_count; i++)
{
if (x[i] <= LP_EPSILON) continue;
log_debug(" x%-5d = %.6f\n", i, x[i]);
@ -404,7 +385,7 @@ int build_edge_map(struct GTSP *gtsp, int *edge_map)
{
int node_count = gtsp->graph->node_count;
int k = node_count;
int k = 0;
for (int i = 0; i < node_count; i++)
{
for (int j = i + 1; j < node_count; j++)
@ -437,7 +418,7 @@ int GTSP_check_solution(struct GTSP *data, double *x)
stack = (struct Node **) malloc(graph->node_count * sizeof(struct Node *));
abort_if(!stack, "could not allocate stack");
for (int i = 0; i < node_count + edge_count; i++)
for (int i = 0; i < edge_count; i++)
{
abort_iff(x[i] < 1.0 - LP_EPSILON && x[i] > LP_EPSILON,
"solution is not integral: x%d = %.4lf", i, x[i]);
@ -453,13 +434,13 @@ int GTSP_check_solution(struct GTSP *data, double *x)
cluster_mark[i] = 0;
int initial;
for (initial = 0; initial < node_count; initial++)
for (initial = 0; initial < edge_count; initial++)
if (x[initial] > 1.0 - LP_EPSILON) break;
abort_if(initial == node_count, "no initial node");
abort_if(initial == edge_count, "no initial node");
stack[stack_top++] = &graph->nodes[initial];
graph->nodes[initial].mark = 1;
stack[stack_top++] = graph->edges[initial].from;
graph->edges[initial].from->mark = 1;
while (stack_top > 0)
{
@ -472,7 +453,7 @@ int GTSP_check_solution(struct GTSP *data, double *x)
struct Node *neighbor = adj->neighbor;
if (neighbor->mark) continue;
if (x[node_count + adj->edge->index] < LP_EPSILON) continue;
if (x[adj->edge->index] < LP_EPSILON) continue;
stack[stack_top++] = neighbor;
neighbor->mark = 1;
@ -504,7 +485,7 @@ int GTSP_solution_found(struct BNC *bnc, struct GTSP *data, double *x)
struct Tour* tour;
tour = (struct Tour*) malloc(data->cluster_count*sizeof(struct Tour));
sprintf(filename, "tmp/gtsp-m%d-n%d-s%d.out", data->cluster_count,
data->graph->node_count, SEED);
@ -517,7 +498,7 @@ int GTSP_solution_found(struct BNC *bnc, struct GTSP *data, double *x)
rval = large_neighborhood_search(tour, data, &tour_cost);
abort_if(rval, "large_neighborhood_search failed");
if(tour_cost + LP_EPSILON < *best_val){
log_info("Local search improve the integral solution\n");
log_info(" obj val = %f\n",*best_val );
@ -684,8 +665,8 @@ int GTSP_main(int argc, char **argv)
rval = GTSP_write_problem(&data, "gtsp.in");
#endif
int init_val;
int init_val = 0;
//
initial_x = (double *) malloc(
(data.graph->node_count + data.graph->edge_count) * sizeof(double));
abort_if(!initial_x, "could not allocate initial_x");
@ -729,7 +710,7 @@ int GTSP_main(int argc, char **argv)
for (int i = 0; i < data.graph->edge_count; i++)
{
struct Edge *e = &data.graph->edges[i];
opt_val += OPTIMAL_X[i + input_node_count] * e->weight;
opt_val += OPTIMAL_X[i] * e->weight;
}
log_info(" opt = %.2lf\n", opt_val);
@ -839,7 +820,7 @@ int build_x_from_tour(struct GTSP *data, struct Tour *tour, double *x)
rval = build_edge_map(data, edge_map);
abort_if(rval, "build_edge_map failed");
for (int i = 0; i < node_count + edge_count; i++)
for (int i = 0; i < edge_count; i++)
x[i] = 0.0;
int next_vertex = tour[0].next;
@ -851,8 +832,6 @@ int build_x_from_tour(struct GTSP *data, struct Tour *tour, double *x)
current_vertex = tour[next_vertex].vertex;
next_vertex = tour[next_vertex].next;
x[from] = 1.0;
x[to] = 1.0;
x[edge_map[from * node_count + to]] = 1.0;
}
@ -1248,6 +1227,7 @@ int build_tour_from_x(struct GTSP *data, struct Tour *tour, double *x)
struct Graph *graph = data->graph;
const int node_count = graph->node_count;
const int edge_count = graph->edge_count;
cluster_mark = (int *) malloc(data->cluster_count * sizeof(int));
abort_if(!cluster_mark, "could not allocate cluster_mark");
@ -1262,10 +1242,12 @@ int build_tour_from_x(struct GTSP *data, struct Tour *tour, double *x)
cluster_mark[i] = 0;
int initial;
for (initial = 0; initial < node_count; initial++)
for (initial = 0; initial < edge_count; initial++)
if (x[initial] > 1.0 - LP_EPSILON) break;
abort_if(initial == node_count, "no initial node");
initial = graph->edges[initial].from->index;
abort_if(initial == edge_count, "no initial node");
stack[stack_top++] = &graph->nodes[initial];
graph->nodes[initial].mark = 1;
@ -1284,7 +1266,7 @@ int build_tour_from_x(struct GTSP *data, struct Tour *tour, double *x)
struct Node *neighbor = adj->neighbor;
if (neighbor->mark) continue;
if (x[node_count + adj->edge->index] < LP_EPSILON) continue;
if (x[adj->edge->index] < LP_EPSILON) continue;
stack[stack_top++] = neighbor;
tour[next_vertex].vertex = neighbor->index;

2
src/params.h
Unescape View File

@ -10,7 +10,7 @@
#define ENABLE_COMB_INEQUALITIES
#define MAX_TOTAL_TIME 3600
#define MAX_TOTAL_TIME 999999
//#define ALLOW_FRACTIONAL_SOLUTIONS