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

Split gtsp.c

Browse Source
This commit is contained in:
Alinson S. Xavier
2015-03-26 18:21:28 -04:00
parent 461aa31475
commit 6e4d69e289
8 changed files with 599 additions and 519 deletions
Download Patch File Download Diff File Expand all files Collapse all files

5
src/flow.c
View File

@@ -54,6 +54,8 @@ int flow_mark_reachable_nodes(
return rval;
}
extern double FLOW_CPU_TIME;
int flow_find_max_flow(
const struct Graph *digraph,
const double *capacities,
@@ -65,6 +67,7 @@ int flow_find_max_flow(
int rval = 0;
FLOW_MAX_FLOW_COUNT++;
double initial_time = get_current_time();
for (int i = 0; i < digraph->node_count; i++)
digraph->nodes[i].mark = 0;
@@ -154,6 +157,8 @@ int flow_find_max_flow(
rval = flow_mark_reachable_nodes(digraph, residual_caps, from);
abort_if(rval, "flow_mark_reachable_nodes failed");
FLOW_CPU_TIME += get_current_time() - initial_time;
CLEANUP:
if (path_edges) free(path_edges);
if (residual_caps) free(residual_caps);

513
src/gtsp-subtour.c Normal file
View File

@@ -0,0 +1,513 @@
#include "gtsp-subtour.h"
#include <assert.h>
#include <float.h>
#include "util.h"
#include "flow.h"
extern double FLOW_CPU_TIME;
int static build_flow_digraph(
struct GTSP *data, double *x, struct Graph *digraph, double *capacities)
{
int rval = 0;
int *digraph_edges = 0;
int node_count = data->graph->node_count;
struct Graph *graph = data->graph;
int digraph_node_count = node_count + data->cluster_count + 1;
int digraph_edge_count = 4 * graph->edge_count + 2 * graph->node_count +
2 * data->cluster_count;
digraph_edges = (int *) malloc(2 * digraph_edge_count * sizeof(int));
abort_if(!digraph_edges, "could not allocate digraph_edges");
// Create four directed edges for each edge of the original graph
int ke = 0;
int kc = 0;
for (int i = 0; i < graph->edge_count; i++)
{
if (x[node_count + i] < LP_EPSILON) continue;
struct Edge *e = &graph->edges[i];
int from = e->from->index;
int to = e->to->index;
digraph_edges[ke++] = from;
digraph_edges[ke++] = to;
capacities[kc++] = x[node_count + i];
digraph_edges[ke++] = to;
digraph_edges[ke++] = from;
capacities[kc++] = 0;
digraph_edges[ke++] = to;
digraph_edges[ke++] = from;
capacities[kc++] = x[node_count + i];
digraph_edges[ke++] = from;
digraph_edges[ke++] = to;
capacities[kc++] = 0;
}
// Create an extra node for each cluster and connect it to the vertices
// of the cluster through some edge with very high capacity
for (int i = 0; i < node_count; i++)
{
struct Node *n = &graph->nodes[i];
int cl = data->clusters[n->index];
digraph_edges[ke++] = n->index;
digraph_edges[ke++] = node_count + cl;
capacities[kc++] = 0;
digraph_edges[ke++] = node_count + cl;
digraph_edges[ke++] = n->index;
capacities[kc++] = 0;
}
// Create an extra root node and connect it to each cluster node through
// some edge with zero capacity
for (int i = 0; i < data->cluster_count; i++)
{
digraph_edges[ke++] = node_count + i;
digraph_edges[ke++] = node_count + data->cluster_count;
capacities[kc++] = 0;
digraph_edges[ke++] = node_count + data->cluster_count;
digraph_edges[ke++] = node_count + i;
capacities[kc++] = 0;
}
assert(ke <= 2 * digraph_edge_count);
assert(kc <= digraph_edge_count);
digraph_edge_count = kc;
rval = graph_build(digraph_node_count, kc, digraph_edges, 1, digraph);
abort_if(rval, "graph_build failed");
for (int i = 0; i < digraph_edge_count; i += 2)
{
digraph->edges[i].reverse = &digraph->edges[i + 1];
digraph->edges[i + 1].reverse = &digraph->edges[i];
}
CLEANUP:
if (digraph_edges) free(digraph_edges);
return rval;
}
int static 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 rval = 0;
char sense = 'G';
double rhs = 2.0 - 2.0 * type;
int newnz = cut_edges_count + type;
int rmatbeg = 0;
int *rmatind = 0;
double *rmatval = 0;
rmatind = (int *) malloc(newnz * sizeof(int));
abort_if(!rmatind, "could not allocate rmatind");
rmatval = (double *) malloc(newnz * sizeof(double));
abort_if(!rmatval, "could not allocate rmatval");
for (int i = 0; i < cut_edges_count; i++)
{
rmatind[i] = cut_edges[i]->index + graph->node_count;
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]);
log_verbose(" %c %.2lf\n", sense, rhs);
if (OPTIMAL_X)
{
double sum = 0;
for (int i = 0; i < newnz; i++)
sum += rmatval[i] * OPTIMAL_X[rmatind[i]];
abort_if(sum <= rhs - LP_EPSILON, "cannot add invalid cut");
}
rval = LP_add_rows(lp, 1, newnz, &rhs, &sense, &rmatbeg, rmatind, rmatval);
abort_if(rval, "LP_add_rows failed");
CLEANUP:
if (rmatval) free(rmatval);
if (rmatind) free(rmatind);
return rval;
}
int find_exact_subtour_cuts(
struct LP *lp, struct GTSP *data, int *total_added_cuts)
{
int rval = 0;
double *x = 0;
double *capacities = 0;
int added_cuts_count = 0;
struct Graph *graph = data->graph;
int num_cols = LP_get_num_cols(lp);
x = (double *) malloc(num_cols * sizeof(double));
abort_if(!x, "could not allocate x");
rval = LP_get_x(lp, x);
abort_if(rval, "LP_get_x failed");
#if LOG_LEVEL >= LOG_LEVEL_DEBUG
rval = GTSP_write_solution(data, "gtsp-frac.out", x);
abort_if(rval, "GTSP_write_solution failed");
#endif
struct Graph digraph;
graph_init(&digraph);
int digraph_edge_count = 4 * graph->edge_count + 2 * graph->node_count +
2 * data->cluster_count;
capacities = (double *) malloc(digraph_edge_count * sizeof(double));
abort_if(!capacities, "could not allocate capacities");
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, &added_cuts_count);
abort_if(rval, "find_exact_subtour_cuts_cluster_to_cluster failed");
log_debug("Added %d cluster-to-cluster subtour cuts\n", added_cuts_count);
(*total_added_cuts) += added_cuts_count;
if (added_cuts_count > 0) goto CLEANUP;
// Constraints (2.2)
rval = find_exact_subtour_cuts_node_to_cluster(lp, data, x, &digraph,
capacities, &added_cuts_count);
abort_if(rval, "find_exact_subtour_cuts_node_to_cluster failed");
log_debug("Added %d node-to-cluster subtour cuts\n", added_cuts_count);
(*total_added_cuts) += added_cuts_count;
if (added_cuts_count > 0) goto CLEANUP;
// Constraints (2.3)
rval = find_exact_subtour_cuts_node_to_node(lp, data, x, &digraph,
capacities, &added_cuts_count);
abort_if(rval, "find_exact_subtour_cuts_node_to_node failed");
log_debug("Added %d node-to-node subtour cuts\n", added_cuts_count);
(*total_added_cuts) += added_cuts_count;
CLEANUP:
graph_free(&digraph);
if (capacities) free(capacities);
if (x) free(x);
return rval;
}
int find_exact_subtour_cuts_node_to_node(
struct LP *lp,
struct GTSP *data,
double *x,
struct Graph *digraph,
double *capacities,
int *added_cuts_count)
{
int rval = 0;
struct Edge **cut_edges = 0;
double *flow = 0;
struct Graph *graph = data->graph;
int *clusters = data->clusters;
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");
(*added_cuts_count)++;
}
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,
int *added_cuts_count)
{
int rval = 0;
int cuts_count = 0;
struct Edge **cut_edges = 0;
double *flow = 0;
struct Graph *graph = data->graph;
int *clusters = data->clusters;
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");
(*added_cuts_count)++;
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,
struct Graph *digraph,
double *capacities,
int *added_cuts_count)
{
int rval = 0;
double *flow = 0;
struct Edge **cut_edges = 0;
int cuts_count = 0;
struct Graph *graph = data->graph;
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");
struct Node *root_node = &digraph->nodes[graph->node_count +
data->cluster_count];
for (int i = 0; i < data->cluster_count; i++)
{
for (int j = i + 1; j < data->cluster_count; j++)
{
struct Node *from = &digraph->nodes[graph->node_count + i];
struct Node *to = &digraph->nodes[graph->node_count + j];
double flow_value;
int cut_edges_count;
activate_cluster_node(capacities, from);
activate_cluster_node(capacities, to);
deactivate_cluster_node(capacities, root_node);
log_verbose("Sending flow from cluster %d to cluster %d\n", i, j);
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, from);
deactivate_cluster_node(capacities, to);
if (flow_value >= 2 - 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, 0, 0, cut_edges, cut_edges_count,
0);
abort_if(rval, "add_subtour_cut failed");
(*added_cuts_count)++;
cuts_count++;
}
// if(cuts_count > 0) break;
}
CLEANUP:
if (cut_edges) free(cut_edges);
if (flow) free(flow);
return rval;
}
void deactivate_cluster_node(double *capacities, struct Node *cluster_node)
{
for (int i = 0; i < cluster_node->degree; i++)
{
struct Adjacency *adj = &cluster_node->adj[i];
struct Edge *e = adj->edge;
capacities[e->index] = 0;
}
}
void activate_cluster_node(double *capacities, struct Node *cluster_node)
{
for (int i = 0; i < cluster_node->degree; i++)
{
struct Adjacency *adj = &cluster_node->adj[i];
struct Edge *e = adj->edge;
capacities[e->index] = 1e10;
capacities[e->reverse->index] = 1e10;
}
}

38
src/gtsp-subtour.h Normal file
View File

@@ -0,0 +1,38 @@
#ifndef _PROJECT_GTSP_SUBTOUR_H_
#define _PROJECT_GTSP_SUBTOUR_H_
#include "gtsp.h"
void deactivate_cluster_node(double *capacities, struct Node *cluster_node);
void activate_cluster_node(double *capacities, struct Node *cluster_node);
int find_exact_subtour_cuts_cluster_to_cluster(
struct LP *lp,
struct GTSP *data,
struct Graph *digraph,
double *capacities,
int *added_cuts_count);
int find_exact_subtour_cuts_node_to_cluster(
struct LP *lp,
struct GTSP *data,
double *x,
struct Graph *digraph,
double *capacities,
int *added_cuts_count);
int find_exact_subtour_cuts_node_to_node(
struct LP *lp,
struct GTSP *data,
double *x,
struct Graph *digraph,
double *capacities,
int *added_cuts_count);
int find_exact_subtour_cuts(
struct LP *lp, struct GTSP *data, int *total_added_cuts);
#endif //_PROJECT_GTSP_SUBTOUR_H_

515
src/gtsp.c
View File

@@ -1,7 +1,6 @@
#include <stdio.h>
#include <stdlib.h>
#include <float.h>
#include <assert.h>
#include <getopt.h>
#include "gtsp.h"
#include "geometry.h"
@@ -9,8 +8,9 @@
#include "flow.h"
#include "branch_and_cut.h"
#include "math.h"
#include "gtsp-subtour.h"
static double *OPTIMAL_X = 0;
double *OPTIMAL_X = 0;
static int get_edge_num(int node_count, int from, int to)
{
@@ -183,507 +183,6 @@ int GTSP_init_lp(struct LP *lp, struct GTSP *data)
return rval;
}
static 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 rval = 0;
char sense = 'G';
double rhs = 2.0 - 2.0 * type;
int newnz = cut_edges_count + type;
int rmatbeg = 0;
int *rmatind = 0;
double *rmatval = 0;
rmatind = (int *) malloc(newnz * sizeof(int));
abort_if(!rmatind, "could not allocate rmatind");
rmatval = (double *) malloc(newnz * sizeof(double));
abort_if(!rmatval, "could not allocate rmatval");
for (int i = 0; i < cut_edges_count; i++)
{
rmatind[i] = cut_edges[i]->index + graph->node_count;
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]);
log_verbose(" %c %.2lf\n", sense, rhs);
if (OPTIMAL_X)
{
double sum = 0;
for (int i = 0; i < newnz; i++)
sum += rmatval[i] * OPTIMAL_X[rmatind[i]];
abort_if(sum <= rhs - LP_EPSILON, "cannot add invalid cut");
}
rval = LP_add_rows(lp, 1, newnz, &rhs, &sense, &rmatbeg, rmatind, rmatval);
abort_if(rval, "LP_add_rows failed");
CLEANUP:
if (rmatval) free(rmatval);
if (rmatind) free(rmatind);
return rval;
}
static int build_flow_digraph(
struct GTSP *data, double *x, struct Graph *digraph, double *capacities)
{
int rval = 0;
int *digraph_edges = 0;
int node_count = data->graph->node_count;
struct Graph *graph = data->graph;
int digraph_node_count = node_count + data->cluster_count + 1;
int digraph_edge_count = 4 * graph->edge_count + 2 * graph->node_count +
2 * data->cluster_count;
digraph_edges = (int *) malloc(2 * digraph_edge_count * sizeof(int));
abort_if(!digraph_edges, "could not allocate digraph_edges");
// Create four directed edges for each edge of the original graph
int ke = 0;
int kc = 0;
for (int i = 0; i < graph->edge_count; i++)
{
// if (x[node_count + i] < LP_EPSILON) continue;
struct Edge *e = &graph->edges[i];
int from = e->from->index;
int to = e->to->index;
digraph_edges[ke++] = from;
digraph_edges[ke++] = to;
capacities[kc++] = x[node_count + i];
digraph_edges[ke++] = to;
digraph_edges[ke++] = from;
capacities[kc++] = 0;
digraph_edges[ke++] = to;
digraph_edges[ke++] = from;
capacities[kc++] = x[node_count + i];
digraph_edges[ke++] = from;
digraph_edges[ke++] = to;
capacities[kc++] = 0;
}
// Create an extra node for each cluster and connect it to the vertices
// of the cluster through some edge with very high capacity
for (int i = 0; i < node_count; i++)
{
struct Node *n = &graph->nodes[i];
int cl = data->clusters[n->index];
digraph_edges[ke++] = n->index;
digraph_edges[ke++] = node_count + cl;
capacities[kc++] = 0;
digraph_edges[ke++] = node_count + cl;
digraph_edges[ke++] = n->index;
capacities[kc++] = 0;
}
// Create an extra root node and connect it to each cluster node through
// some edge with zero capacity
for (int i = 0; i < data->cluster_count; i++)
{
digraph_edges[ke++] = node_count + i;
digraph_edges[ke++] = node_count + data->cluster_count;
capacities[kc++] = 0;
digraph_edges[ke++] = node_count + data->cluster_count;
digraph_edges[ke++] = node_count + i;
capacities[kc++] = 0;
}
assert(ke <= 2 * digraph_edge_count);
assert(kc <= digraph_edge_count);
// digraph_edge_count = kc;
rval = graph_build(digraph_node_count, kc, digraph_edges, 1, digraph);
abort_if(rval, "graph_build failed");
for (int i = 0; i < digraph_edge_count; i += 2)
{
digraph->edges[i].reverse = &digraph->edges[i + 1];
digraph->edges[i + 1].reverse = &digraph->edges[i];
}
CLEANUP:
if (digraph_edges) free(digraph_edges);
return rval;
}
void deactivate_cluster_node(double *capacities, struct Node *cluster_node)
{
for (int i = 0; i < cluster_node->degree; i++)
{
struct Adjacency *adj = &cluster_node->adj[i];
struct Edge *e = adj->edge;
capacities[e->index] = 0;
}
}
void activate_cluster_node(double *capacities, struct Node *cluster_node)
{
for (int i = 0; i < cluster_node->degree; i++)
{
struct Adjacency *adj = &cluster_node->adj[i];
struct Edge *e = adj->edge;
capacities[e->index] = 1e10;
capacities[e->reverse->index] = 1e10;
}
}
int find_exact_subtour_cuts_cluster_to_cluster(
struct LP *lp,
struct GTSP *data,
struct Graph *digraph,
double *capacities,
int *added_cuts_count)
{
int rval = 0;
double *flow = 0;
struct Edge **cut_edges = 0;
int cuts_count = 0;
struct Graph *graph = data->graph;
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");
struct Node *root_node = &digraph->nodes[graph->node_count +
data->cluster_count];
for (int i = 0; i < data->cluster_count; i++)
{
for (int j = i + 1; j < data->cluster_count; j++)
{
struct Node *from = &digraph->nodes[graph->node_count + i];
struct Node *to = &digraph->nodes[graph->node_count + j];
double flow_value;
int cut_edges_count;
activate_cluster_node(capacities, from);
activate_cluster_node(capacities, to);
deactivate_cluster_node(capacities, root_node);
log_verbose("Sending flow from cluster %d to cluster %d\n", i, j);
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, from);
deactivate_cluster_node(capacities, to);
if (flow_value >= 2 - LP_EPSILON) 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, 0, 0, cut_edges, cut_edges_count,
0);
abort_if(rval, "add_subtour_cut failed");
(*added_cuts_count)++;
cuts_count++;
}
}
CLEANUP:
if (cut_edges) free(cut_edges);
if (flow) free(flow);
return rval;
}
int find_exact_subtour_cuts_node_to_cluster(
struct LP *lp,
struct GTSP *data,
double *x,
struct Graph *digraph,
double *capacities,
int *added_cuts_count)
{
int rval = 0;
int cuts_count = 0;
struct Edge **cut_edges = 0;
double *flow = 0;
struct Graph *graph = data->graph;
int *clusters = data->clusters;
cut_edges = (struct Edge **) malloc(
digraph->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 + LP_EPSILON >= 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");
(*added_cuts_count)++;
cuts_count++;
}
}
CLEANUP:
if (cut_edges) free(cut_edges);
if (flow) free(flow);
return rval;
}
int find_exact_subtour_cuts_node_to_node(
struct LP *lp,
struct GTSP *data,
double *x,
struct Graph *digraph,
double *capacities,
int *added_cuts_count)
{
int rval = 0;
struct Edge **cut_edges = 0;
double *flow = 0;
struct Graph *graph = data->graph;
int *clusters = data->clusters;
cut_edges = (struct Edge **) malloc(
digraph->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) - LP_EPSILON)
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");
(*added_cuts_count)++;
}
CLEANUP:
if (flow) free(flow);
if (cut_edges) free(cut_edges);
return rval;
}
int find_exact_subtour_cuts(
struct LP *lp, struct GTSP *data, int *total_added_cuts)
{
int rval = 0;
double *x = 0;
double *capacities = 0;
int added_cuts_count = 0;
struct Graph *graph = data->graph;
int num_cols = LP_get_num_cols(lp);
x = (double *) malloc(num_cols * sizeof(double));
abort_if(!x, "could not allocate x");
rval = LP_get_x(lp, x);
abort_if(rval, "LP_get_x failed");
#if LOG_LEVEL >= LOG_LEVEL_DEBUG
rval = GTSP_write_solution(data, "gtsp-frac.out", x);
abort_if(rval, "GTSP_write_solution failed");
#endif
struct Graph digraph;
graph_init(&digraph);
int digraph_edge_count = 4 * graph->edge_count + 2 * graph->node_count +
2 * data->cluster_count;
capacities = (double *) malloc(digraph_edge_count * sizeof(double));
abort_if(!capacities, "could not allocate capacities");
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, &added_cuts_count);
abort_if(rval, "find_exact_subtour_cuts_cluster_to_cluster failed");
log_debug("Added %d cluster-to-cluster subtour cuts\n", added_cuts_count);
(*total_added_cuts) += added_cuts_count;
if (added_cuts_count > 0) goto CLEANUP;
// Constraints (2.2)
rval = find_exact_subtour_cuts_node_to_cluster(lp, data, x, &digraph,
capacities, &added_cuts_count);
abort_if(rval, "find_exact_subtour_cuts_node_to_cluster failed");
log_debug("Added %d node-to-cluster subtour cuts\n", added_cuts_count);
(*total_added_cuts) += added_cuts_count;
// Constraints (2.3)
rval = find_exact_subtour_cuts_node_to_node(lp, data, x, &digraph,
capacities, &added_cuts_count);
abort_if(rval, "find_exact_subtour_cuts_node_to_node failed");
log_debug("Added %d node-to-node subtour cuts\n", added_cuts_count);
(*total_added_cuts) += added_cuts_count;
CLEANUP:
graph_free(&digraph);
if (capacities) free(capacities);
if (x) free(x);
return rval;
}
int GTSP_add_cutting_planes(struct LP *lp, struct GTSP *data)
{
int rval = 0;
@@ -952,6 +451,11 @@ int GTSP_solution_found(struct GTSP *data, double *x)
return rval;
}
double FLOW_CPU_TIME = 0;
double LP_CPU_TIME = 0;
int LP_OPTIMIZE_COUNT = 0;
int GTSP_main(int argc, char **argv)
{
int rval = 0;
@@ -1026,7 +530,10 @@ int GTSP_main(int argc, char **argv)
log_info(" obj value = %.2lf **\n", bnc.best_obj_val);
log_info("Branch-and-bound nodes: %d\n", BNC_NODE_COUNT);
log_info("Max-flow computations: %d\n", FLOW_MAX_FLOW_COUNT);
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);
CLEANUP:
GTSP_free(&data);

5
src/gtsp.h
View File

@@ -19,6 +19,8 @@ struct GTSP
double *y_coordinates;
};
static const double MIN_CUT_VIOLATION = 0.5;
int GTSP_create_random_problem(
int node_count, int cluster_count, int grid_size, struct GTSP *data);
@@ -36,4 +38,7 @@ int GTSP_write_solution(struct GTSP *data, char *filename, double *x);
int GTSP_main(int argc, char **argv);
extern double *OPTIMAL_X;
extern double FLOW_CPU_TIME;
#endif //_PROJECT_GTSP_H_

13
src/lp.c
View File

@@ -119,25 +119,32 @@ int LP_change_bound(struct LP *lp, int col, char lower_or_upper, double bnd)
return rval;
}
extern int LP_OPTIMIZE_COUNT;
extern double LP_CPU_TIME;
int LP_optimize(struct LP *lp, int *infeasible)
{
LP_OPTIMIZE_COUNT++;
int rval = 0, solstat;
*infeasible = 0;
double current = get_current_time();
rval = CPXdualopt(lp->cplex_env, lp->cplex_lp);
abort_if(rval, "CPXdualopt failed");
LP_CPU_TIME += get_current_time() - current;
solstat = CPXgetstat(lp->cplex_env, lp->cplex_lp);
if (solstat == CPX_STAT_INFEASIBLE)
{
*infeasible = 1;
goto CLEANUP;
}
else
{
abort_if(solstat != CPX_STAT_OPTIMAL
&& solstat != CPX_STAT_OPTIMAL_INFEAS,
"Invalid solution status");
abort_if(solstat != CPX_STAT_OPTIMAL &&
solstat != CPX_STAT_OPTIMAL_INFEAS, "Invalid solution status");
}
CLEANUP:

21
src/tsp.c
View File

@@ -64,8 +64,7 @@ int TSP_init_lp(struct LP *lp, struct TSPData *data)
}
int TSP_find_violated_subtour_elimination_cut(
struct LP *lp,
struct TSPData *data)
struct LP *lp, struct TSPData *data)
{
int ncount = data->node_count;
int edge_count = data->edge_count;
@@ -118,9 +117,14 @@ int TSP_find_violated_subtour_elimination_cut(
int delta_count = 0;
int island_count = 0;
while (!TSP_is_graph_connected(&G, x, &island_count, island_sizes,
island_start, island_nodes))
while (1)
{
rval = graph_find_connected_components(&G, x, &island_count,
island_sizes, island_start, island_nodes);
abort_if(rval, "graph_find_connected_components failed");
if(island_count == 1) break;
log_verbose("Adding %d subtour inequalities...\n", island_count);
for (int i = 0; i < island_count; i++)
{
@@ -184,7 +188,7 @@ int TSP_add_subtour_elimination_cut(struct LP *lp, int delta_length, int *delta)
return rval;
}
int TSP_is_graph_connected(
int graph_find_connected_components(
struct Graph *G,
double *x,
int *island_count,
@@ -216,7 +220,7 @@ int TSP_is_graph_connected(
(*island_count) = current_island;
return (*island_count == 1);
return 0;
}
int TSP_find_closest_neighbor_tour(
@@ -491,8 +495,9 @@ int TSP_main(int argc, char **argv)
bnc.best_obj_val = TSP_find_initial_solution(&data);
bnc.problem_data = (void *) &data;
bnc.problem_init_lp = (int (*)(struct LP *, void *)) TSP_init_lp;
bnc.problem_add_cutting_planes =
(int (*)(struct LP *, void *)) TSP_add_cutting_planes;
bnc.problem_add_cutting_planes = (int (*)(
struct LP *,
void *)) TSP_add_cutting_planes;
rval = BNC_init_lp(&bnc);
abort_if(rval, "BNC_init_lp failed");

2
src/tsp.h
View File

@@ -19,7 +19,7 @@ void TSP_free_data(struct TSPData *data);
int TSP_find_violated_subtour_elimination_cut
(struct LP *lp, struct TSPData *data);
int TSP_is_graph_connected(
int graph_find_connected_components(
struct Graph *G,
double *x,
int *island_count,