1 /*===========================================================================*/
2 /* */
3 /* This file is part of a demonstration application for use with the */
4 /* SYMPHONY Branch, Cut, and Price Library. This application is a solver for */
5 /* Capacitated Network Routing Problems. */
6 /* */
7 /* (c) Copyright 2000-2007 Ted Ralphs. All Rights Reserved. */
8 /* */
9 /* This application was developed by Ted Ralphs (ted@lehigh.edu) */
10 /* */
11 /* This software is licensed under the Eclipse Public License. Please see */
12 /* accompanying file for terms. */
13 /* */
14 /*===========================================================================*/
15
16 /* system include files */
17 #include <stdlib.h>
18 #include <math.h>
19 #include <string.h> /* memset() is defined here in LINUX */
20
21 /* SYMPHONY include files */
22 #include "sym_constants.h"
23 #include "sym_macros.h"
24 #include "sym_messages.h"
25 #include "sym_proccomm.h"
26 #include "sym_cg.h"
27
28 /* CNRP include files */
29 #include "cnrp_cg.h"
30 #include "cnrp_macros.h"
31 #include "network.h"
32
33 /*__BEGIN_EXPERIMENTAL_SECTION__*/
34 #if 0
35 #include "util.h"
36 extern CCrandstate rand_state;
37 #endif
38 /*___END_EXPERIMENTAL_SECTION___*/
39
40 /*===========================================================================*/
41
42 #define SEND_DIR_SUBTOUR_CONSTRAINT(num_nodes, total_demand) \
43 new_cut->type = (num_nodes < vertnum/2 ? \
44 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS); \
45 new_cut->rhs = (new_cut->type == SUBTOUR_ELIM_SIDE ? \
46 RHS(num_nodes, total_demand, capacity) : \
47 BINS(total_demand, capacity)); \
48 cuts_found += cg_send_cut(new_cut, num_cuts, alloc_cuts, cuts); \
49
50 #define SEND_SUBTOUR_CONSTRAINT(num_nodes, total_demand) \
51 if (mult - 1){ \
52 new_cut->type = (num_nodes < vertnum/2 ? \
53 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS); \
54 new_cut->rhs = (new_cut->type == SUBTOUR_ELIM_SIDE ? \
55 RHS(num_nodes, total_demand, capacity) : \
56 mult*BINS(total_demand, capacity)); \
57 cuts_found += cg_send_cut(new_cut, num_cuts, alloc_cuts, cuts); \
58 }else{ \
59 new_cut->type = SUBTOUR_ELIM_ACROSS; \
60 new_cut->rhs = mult*BINS(total_demand, capacity); \
61 cuts_found += cg_send_cut(new_cut, num_cuts, alloc_cuts, cuts); \
62 } \
63
64 /*===========================================================================*/
65
66 /*===========================================================================*\
67 * This file implements the greedy shrinking algorithm of Augerat, et al.
68 * The original implementation was done by Leonid Kopman.
69 \*===========================================================================*/
70
reduce_graph(network * n,double etol,double * demand,double capacity,int mult,cut_data * new_cut,int * num_cuts,int * alloc_cuts,cut_data *** cuts)71 int reduce_graph(network *n, double etol, double *demand, double capacity,
72 int mult, cut_data *new_cut, int *num_cuts, int *alloc_cuts,
73 cut_data ***cuts)
74 {
75 elist *e1, *e2, *e3;
76 edge *cur_edge;
77 int v1, v2, deg, count, i, k, vertnum = n->vertnum;
78 edge *edges = n->edges;
79 int num_edges = n->edgenum, cuts_found = 0;
80 int edges_deleted = 0;
81 vertex *verts = n->verts;
82 vertex *v2_pt, *third_node;
83 char *coef;
84
85 new_cut->size = (vertnum >> DELETE_POWER) + 1;
86 new_cut->coef = coef = (char *) (calloc(new_cut->size, CSIZE));
87
88 while(TRUE){
89 num_edges = n->edgenum;
90 edges_deleted = 0;
91 for (i = 0; i < num_edges; i++){
92 cur_edge = edges + i;
93 if (cur_edge->weight >= 1 - etol && cur_edge->v0 &&
94 cur_edge->v1 && !(cur_edge->deleted)){
95 cur_edge->deleted = TRUE;
96 n->edgenum--;
97 edges_deleted++;
98 v1 = (verts[cur_edge->v0].degree ==
99 MIN(verts[cur_edge->v0].degree, verts[cur_edge->v1].degree))?
100 cur_edge->v0 : cur_edge->v1;
101 v2 = (v1 == cur_edge->v0) ? cur_edge->v1 : cur_edge->v0;
102 if (cur_edge->weight + verts[v1].weight + verts[v2].weight >
103 RHS(verts[v1].orig_node_list_size +
104 verts[v2].orig_node_list_size + 2,
105 demand[v1]+demand[v2], capacity)){
106 memset(coef, 0, new_cut->size*CSIZE);
107 for (k = 0; k < verts[v1].orig_node_list_size; k++)
108 (coef[(verts[v1].orig_node_list)[k] >>
109 DELETE_POWER]) |=
110 (1 << ((verts[v1].orig_node_list)[k] &
111 DELETE_AND));
112 (coef[v1 >> DELETE_POWER]) |=
113 (1 << (v1 & DELETE_AND));
114 for (k = 0; k < verts[v2].orig_node_list_size; k++)
115 (coef[(verts[v2].orig_node_list)[k] >>
116 DELETE_POWER]) |=
117 (1 << ((verts[v2].orig_node_list)[k] &
118 DELETE_AND));
119 (coef[v2 >> DELETE_POWER]) |=
120 (1 << (v2 & DELETE_AND));
121 #if 0
122 new_cut->type = SUBTOUR_ELIM_SIDE;
123 new_cut->rhs = RHS(verts[v1].orig_node_list_size +
124 verts[v2].orig_node_list_size + 2,
125 demand[v1]+demand[v2], capacity);
126 cuts_found += cg_send_cut(new_cut);
127 #endif
128 #ifdef DIRECTED_X_VARS
129 SEND_DIR_SUBTOUR_CONSTRAINT(verts[v1].orig_node_list_size +
130 verts[v2].orig_node_list_size + 2,
131 demand[v1]+demand[v2]);
132 #else
133 SEND_SUBTOUR_CONSTRAINT(verts[v1].orig_node_list_size +
134 verts[v2].orig_node_list_size + 2,
135 demand[v1]+demand[v2]);
136 #endif
137 }
138 verts[v1].deleted = TRUE;
139 demand[v2] += demand[v1];
140 demand[v1] = 0;
141 verts[v2].weight += cur_edge->weight + verts[v1].weight;
142 v2_pt = verts + v2;
143 v2_pt->degree--;
144 if (v2_pt->first->other_end == v1){
145 v2_pt->first = v2_pt->first->next_edge;
146 }else{
147 for (e3 = v2_pt->first; e3 && e3->next_edge; e3 = e3->next_edge)
148 if (e3->next_edge->other_end == v1){
149 e3->next_edge = e3->next_edge->next_edge;
150 if (e3->next_edge == NULL) v2_pt->last = e3;
151 break;
152 }
153 }
154
155 if (!(v2_pt->orig_node_list_size))
156 v2_pt->orig_node_list = (int *) malloc(vertnum*ISIZE);
157 (v2_pt->orig_node_list)[(v2_pt->orig_node_list_size)++] = v1;
158
159 for (k = 0; k < verts[v1].orig_node_list_size; k++){
160 (v2_pt->orig_node_list)[(v2_pt->orig_node_list_size)++] =
161 (verts[v1].orig_node_list)[k];
162 }
163 deg = verts[v1].degree;
164
165 for (e1 = verts[v1].first, count=0; e1 && (count < deg); count++ ){
166 third_node = e1->other;
167 if (third_node->orignodenum == v2){
168 e1 = e1->next_edge;
169 continue;
170 }
171 for (e2 = v2_pt->first; e2; e2 = e2->next_edge){
172 if (e2->other_end == e1->other_end ){
173 e2->data->weight += e1->data->weight;
174 e1->data->deleted = TRUE;
175 edges_deleted++;
176 (third_node->degree)--;
177 if (third_node->first->other_end == v1){
178 third_node->first=third_node->first->next_edge;
179 }else{
180 for (e3 = third_node->first; e3 && e3->next_edge;
181 e3 = e3->next_edge)
182 if (e3->next_edge->other_end == v1){
183 e3->next_edge = e3->next_edge->next_edge;
184 if (e3->next_edge == NULL) third_node->last = e3;
185 break;
186 }
187 }
188 break;
189 }
190 }
191 if (e2){
192 e1 = e1->next_edge;
193 continue;
194 }
195 /* ok, so e1->other_node is not incident to v2 */
196 for (e3 = third_node->first; e3 ; e3 = e3->next_edge){
197 if (e3->other_end == v1){
198 e3->other = v2_pt;
199 e3->other_end = v2;
200 e3->data->v0 = MIN(v2, third_node->orignodenum);
201 e3->data->v1 = MAX(v2, third_node->orignodenum);
202 break;
203 }
204 }
205 v2_pt->last->next_edge = e1;
206 v2_pt->last = e1;
207 v2_pt->degree++;
208 e1=e1->next_edge;
209 v2_pt->last->next_edge = NULL;
210 }
211 }
212 }
213 if (!edges_deleted) break;
214 }
215 FREE(new_cut->coef);
216 return(cuts_found);
217 }
218
219 /*===========================================================================*/
220
greedy_shrinking1(network * n,double capacity,double etol,int max_num_cuts,cut_data * new_cut,int * compnodes,int * compmembers,int compnum,char * in_set,double * cut_val,int * ref,char * cut_list,double * demand,int mult,int * num_cuts,int * alloc_cuts,cut_data *** cuts)221 int greedy_shrinking1(network *n, double capacity, double etol,
222 int max_num_cuts, cut_data *new_cut,
223 int *compnodes, int *compmembers, int compnum,
224 char *in_set, double *cut_val, int *ref, char *cut_list,
225 double *demand, int mult, int *num_cuts, int *alloc_cuts,
226 cut_data ***cuts)
227 {
228 double set_weight, set_demand;
229 vertex *verts = n->verts;
230 elist *e;
231 int cuts_found = 0, i, j, k;
232 char *pt, *cutpt;
233 int *ipt;
234 double *dpt;
235 int vertnum = n->vertnum;
236
237 int max_vert = 0, set_size, begin = 1, cur_comp, end = 1, other_end;
238 char *coef;
239 double maxval, weight;
240 vertex *cur_nodept;
241
242 new_cut->size = (vertnum >> DELETE_POWER) + 1;
243 new_cut->coef = coef = (char *) (calloc(new_cut->size, sizeof(char)));
244 memset(cut_list, 0, new_cut->size * (max_num_cuts + 1));
245
246 *in_set = 0;
247
248 for (i = 1; i < vertnum; i++){
249 if (verts[compmembers[i]].deleted) compmembers[i] = 0;
250 ref[compmembers[i]] = i;
251 }
252 *ref = 0;
253 /* ref is a reference array for compmembers: gives a place
254 in which a vertex is listed in compmembers */
255
256 for (cur_comp = 1; cur_comp <= compnum;
257 begin += compnodes[cur_comp], cur_comp++) /* for every component */
258 for (i = begin, end = begin + compnodes[cur_comp]; i < end; i++){
259 if (compmembers[i] == 0) continue;
260 /* for every node as a starting one */
261 /* initialize the data structures */
262 memset(in_set + begin, 0, compnodes[cur_comp] * sizeof(char));
263 memset(cut_val + begin, 0, compnodes[cur_comp] * sizeof(double));
264 in_set[i] = 1;
265 set_size = 1 + verts[compmembers[i]].orig_node_list_size;
266 set_weight = verts[compmembers[i]].weight;
267 for (e = verts[compmembers[i]].first; e; e = e->next_edge){
268 if (e->other_end)
269 cut_val[ref[e->other_end]] = e->data->weight;
270 }
271 set_demand = demand[compmembers[i]];
272
273 while(TRUE){
274 if (set_weight > RHS(set_size, set_demand, capacity) + etol &&
275 set_size > 2){
276 memset(coef, 0, new_cut->size*sizeof(char));
277 for (j = begin, ipt = compmembers + begin; j < end; j++, ipt++){
278 if (in_set[j]){
279 cur_nodept = verts + (*ipt);
280 if (cur_nodept->orig_node_list_size)
281 for (k = 0; k < cur_nodept->orig_node_list_size; k++)
282 (coef[(cur_nodept->orig_node_list)[k] >>
283 DELETE_POWER]) |=
284 (1 << ((cur_nodept->orig_node_list)[k] &
285 DELETE_AND));
286 (coef[(*ipt) >> DELETE_POWER]) |=
287 (1 << ((*ipt) & DELETE_AND));
288 }
289 }
290 for (k = 0, cutpt = cut_list; k < cuts_found; k++,
291 cutpt += new_cut->size)
292 if (!memcmp(coef, cutpt, new_cut->size*sizeof(char)))
293 break;/* same cuts */
294 if (k >= cuts_found){
295 #if 0
296 new_cut->type = SUBTOUR_ELIM_SIDE;
297 new_cut->rhs = RHS(set_size, set_demand, capacity);
298 cuts_found += cg_send_cut(new_cut);
299 #endif
300 #ifdef DIRECTED_X_VARS
301 SEND_DIR_SUBTOUR_CONSTRAINT(set_size, set_demand);
302 #else
303 SEND_SUBTOUR_CONSTRAINT(set_size, set_demand);
304 #endif
305 memcpy(cutpt, coef, new_cut->size);
306 }
307 if (cuts_found > max_num_cuts){
308 FREE(new_cut->coef);
309 return(cuts_found);
310 }
311 }
312 for (maxval = -1, pt = in_set + begin, dpt = cut_val + begin,
313 j = begin; j < end; pt++, dpt++, j++){
314 if (!(*pt) && *dpt > maxval){
315 maxval = cut_val[j];
316 max_vert = j;
317 }
318 }
319 if (maxval > 0){ /* add the vertex to the set */
320 in_set[max_vert]=1;
321 set_size += 1+ verts[compmembers[max_vert]].orig_node_list_size;
322 set_demand += demand[compmembers[max_vert]];
323 set_weight += verts[compmembers[max_vert]].weight;
324 cut_val[max_vert] = 0;
325 for (e=verts[compmembers[max_vert]].first; e; e = e->next_edge){
326 other_end = ref[e->other_end];
327 weight = e->data->weight;
328 set_weight += (in_set[other_end]) ? weight : -weight;
329 cut_val[other_end] += (in_set[other_end]) ? 0 : weight;
330 }
331 }
332 else{ /* can't add anything to the set */
333 break;
334 }
335 }
336 }
337 FREE(new_cut->coef);
338 return(cuts_found);
339 }
340
341 /*===========================================================================*/
342
343 #if defined(ADD_FLOW_VARS) && defined(DIRECTED_X_VARS)
greedy_shrinking1_dicut(network * n,double capacity,double etol,int max_num_cuts,cut_data * new_cut,int * compnodes,int * compmembers,int compnum,char * in_set,double * cut_val,int * ref,char * cut_list,double * demand,int mult,int * num_cuts,int * alloc_cuts,cut_data *** cuts)344 int greedy_shrinking1_dicut(network *n, double capacity, double etol,
345 int max_num_cuts, cut_data *new_cut,
346 int *compnodes, int *compmembers, int compnum,
347 char *in_set, double *cut_val, int *ref,
348 char *cut_list, double *demand, int mult,
349 int *num_cuts, int *alloc_cuts,
350 cut_data ***cuts)
351 {
352 double set_demand, set_cut_val, tmp_cut_val;
353 vertex *verts = n->verts;
354 elist *e;
355 int cuts_found = 0, i, j;
356 char *pt;
357 int vertnum = n->vertnum;
358
359 int min_vert = 0, set_size, begin = 1, cur_comp, end = 1;
360 char *coef;
361 double min_val;
362 int max_size, numarcs, *arcs;
363
364 max_size = DSIZE + ISIZE + (vertnum >> DELETE_POWER) + 1
365 + 2*vertnum*vertnum*ISIZE;
366 new_cut->coef = (char *) calloc(max_size, sizeof(char));
367 coef = new_cut->coef + DSIZE + ISIZE;
368 arcs = (int *) (coef + (vertnum >> DELETE_POWER) + 1);
369
370 *in_set = 0;
371
372 for (i = 1; i < vertnum; i++){
373 if (verts[compmembers[i]].deleted) compmembers[i] = 0;
374 ref[compmembers[i]] = i;
375 }
376 *ref = 0;
377 /* ref is a reference array for compmembers: gives a place
378 in which a vertex is listed in compmembers */
379
380 for (cur_comp = 1; cur_comp <= compnum; begin += compnodes[cur_comp],
381 cur_comp++){ /* for every component */
382 for (i = begin, end = begin + compnodes[cur_comp]; i < end; i++){
383 if (compmembers[i] == 0) continue;
384 /* for every node as a starting one */
385 /* initialize the data structures */
386 memset(in_set + begin, 0, compnodes[cur_comp] * sizeof(char));
387 in_set[i] = 1;
388 set_size = 1 + verts[compmembers[i]].orig_node_list_size;
389 set_demand = demand[compmembers[i]];
390 set_cut_val = 0;
391 for (e = verts[compmembers[i]].first; e; e = e->next_edge){
392 if (e->other_end < compmembers[i]){
393 set_cut_val += MIN(e->data->flow1,
394 MIN(set_demand, capacity)*e->data->weight1);
395 }else{
396 set_cut_val += MIN(e->data->flow2,
397 MIN(set_demand, capacity)*e->data->weight2);
398 }
399 }
400
401 while(TRUE){
402 if (set_cut_val + etol < set_demand){
403 memset(coef, 0, ((vertnum >> DELETE_POWER) + 1)*sizeof(char));
404 numarcs = 0;
405 for (j = begin; j < end; j++){
406 if (in_set[j]){
407 (coef[(compmembers[j]) >> DELETE_POWER]) |=
408 (1 << ((compmembers[j]) & DELETE_AND));
409 for (e = verts[compmembers[j]].first; e; e=e->next_edge){
410 if (e->other_end < compmembers[j]){
411 if (!in_set[ref[e->other_end]] &&
412 set_demand*e->data->weight1 >
413 e->data->flow1+etol){
414 arcs[numarcs << 1] = e->other_end;
415 arcs[(numarcs << 1) + 1] = compmembers[j];
416 numarcs++;
417 }
418 }else{
419 if (!in_set[ref[e->other_end]] &&
420 set_demand*e->data->weight2 >
421 e->data->flow2+etol){
422 arcs[numarcs << 1] = e->other_end;
423 arcs[(numarcs << 1) + 1] = compmembers[j];
424 numarcs++;
425 }
426 }
427 }
428 }
429 }
430 ((double *)(new_cut->coef))[0] = set_demand;
431 ((int *)(new_cut->coef + DSIZE))[0] = numarcs;
432 new_cut->size = DSIZE + ISIZE + (vertnum >> DELETE_POWER)
433 + 1 + 2 * numarcs * ISIZE;
434 new_cut->type = MIXED_DICUT;
435 new_cut->rhs = set_demand;
436 new_cut->name = CUT__SEND_TO_CP;
437 cuts_found += cg_send_cut(new_cut, num_cuts, alloc_cuts, cuts);
438 if (cuts_found > max_num_cuts){
439 FREE(new_cut->coef);
440 return(cuts_found);
441 }
442 }
443 for (min_val = 0, pt = in_set + begin, j = begin; j < end;
444 pt++, j++){
445 if (!in_set[j]){
446 for (tmp_cut_val = 0, e = verts[compmembers[j]].first;
447 e; e = e->next_edge){
448 if (in_set[ref[e->other_end]]){
449 if (compmembers[j] < e->other_end){
450 tmp_cut_val -=
451 MIN(e->data->flow1,
452 MIN(set_demand, capacity)*e->data->weight1);
453 }else{
454 tmp_cut_val -=
455 MIN(e->data->flow2,
456 MIN(set_demand, capacity)*e->data->weight2);
457 }
458 }else{
459 if (compmembers[j] < e->other_end){
460 tmp_cut_val +=
461 MIN(e->data->flow2,
462 MIN(set_demand, capacity)*e->data->weight2);
463 }else{
464 tmp_cut_val +=
465 MIN(e->data->flow1,
466 MIN(set_demand, capacity)*e->data->weight1);
467 }
468 }
469 }
470 if (tmp_cut_val < min_val - etol){
471 min_vert = j;
472 min_val = tmp_cut_val;
473 }
474 }
475 }
476 if (min_val < 0){ /* add the vertex to the set */
477 in_set[min_vert] = 1;
478 set_size += 1 +
479 verts[compmembers[min_vert]].orig_node_list_size;
480 set_demand += demand[compmembers[min_vert]];
481 set_cut_val += min_val;
482 }
483 else{ /* can't add anything to the set */
484 break;
485 }
486 }
487 }
488 }
489 FREE(new_cut->coef);
490 return(cuts_found);
491 }
492 #endif
493
494 /*===========================================================================*/
495
greedy_shrinking6(network * n,double capacity,double etol,cut_data * new_cut,int * compnodes,int * compmembers,int compnum,char * in_set,double * cut_val,int * ref,char * cut_list,int max_num_cuts,double * demand,int trial_num,double prob,int mult,int * num_cuts,int * alloc_cuts,cut_data *** cuts)496 int greedy_shrinking6(network *n, double capacity, double etol,
497 cut_data *new_cut, int *compnodes,
498 int *compmembers, int compnum,char *in_set,
499 double *cut_val, int *ref, char *cut_list,
500 int max_num_cuts, double *demand, int trial_num,
501 double prob, int mult, int *num_cuts, int *alloc_cuts,
502 cut_data ***cuts)
503 {
504 double set_weight, set_demand;
505 vertex *verts = n->verts;
506 elist *e;
507 int i, j, k, cuts_found = 0;
508 char *pt, *cutpt;
509 double *dpt;
510 int vertnum = n->vertnum;
511
512 int max_vert = 0, set_size, begin = 1, cur_comp, end = 1, num_trials;
513 char *coef;
514 double maxval;
515 double denominator=pow(2.0,31.0)-1.0;
516 double r, q;
517
518 int other_end;
519 double weight;
520 int *ipt;
521 vertex *cur_nodept;
522
523 new_cut->size = (vertnum >> DELETE_POWER) + 1;
524 new_cut->coef =coef= (char *) (calloc(new_cut->size,sizeof(char)));
525 memset(cut_list, 0, new_cut->size * (max_num_cuts +1));
526
527
528 *in_set=0;
529
530 for(i = 1; i < vertnum; i++){
531 if (verts[compmembers[i]].deleted) compmembers[i] = 0;
532 ref[compmembers[i]] = i;
533 }
534 *ref = 0;
535
536 /* ref is a reference array for compmembers: gives a place
537 in which a vertex is listed in compmembers */
538
539 for (cur_comp = 1; cur_comp <= compnum; begin += compnodes[cur_comp],
540 cur_comp++){
541 /* for every component */
542 if (compnodes[cur_comp] <= 7) continue;
543
544 for (num_trials = 0; num_trials < trial_num * compnodes[cur_comp];
545 num_trials++){
546 end = begin + compnodes[cur_comp];
547 /*initialize the data structures */
548 memset(in_set + begin, 0, compnodes[cur_comp] * sizeof(char));
549 memset(cut_val+ begin, 0, compnodes[cur_comp] * sizeof(double));
550 set_size = 0;
551 set_demand = 0;
552 set_weight = 0;
553 for (i = begin; i < end; i++ ){
554 if (compmembers[i] == 0) continue;
555 /*__BEGIN_EXPERIMENTAL_SECTION__*/
556 #if 0
557 r = CCutil_lprand(&rand_state)/CC_PRANDMAX;
558 #endif
559 /*___END_EXPERIMENTAL_SECTION___*/
560 r = RANDOM()/denominator;
561 q = (prob/compnodes[cur_comp]);
562 if (r < q){
563 in_set[i] = 1;
564 set_size += 1 + verts[compmembers[i]].orig_node_list_size;
565 set_demand += demand[compmembers[i]];
566 set_weight += verts[compmembers[i]].weight;
567 for (e = verts[compmembers[i]].first; e; e = e-> next_edge){
568 other_end = ref[e->other_end];
569 weight = e->data->weight;
570 set_weight += (in_set[other_end]) ? weight : -weight;
571 cut_val[other_end] += (in_set[other_end]) ? 0 : weight;
572 }
573 }
574 }
575 while(set_size){
576 if (set_weight > RHS(set_size, set_demand, capacity) + etol &&
577 set_size > 2){
578 memset(coef, 0, new_cut->size*sizeof(char));
579 for (j = begin, ipt = compmembers + begin; j < end; j++, ipt++){
580 if (in_set[j]){
581 cur_nodept = verts + (*ipt);
582 if (cur_nodept->orig_node_list_size)
583 for (k = 0; k < cur_nodept->orig_node_list_size; k++)
584 (coef[(cur_nodept->orig_node_list)[k] >>
585 DELETE_POWER]) |=
586 (1 << ((cur_nodept->orig_node_list)[k] &
587 DELETE_AND));
588 (coef[(*ipt) >> DELETE_POWER]) |= (1 << ((*ipt) &
589 DELETE_AND));
590 }
591 }
592 for (k = 0, cutpt = cut_list; k < cuts_found; k++,
593 cutpt += new_cut->size)
594 if (!memcmp(coef, cutpt, new_cut->size*sizeof(char))) break;
595 if ( k >= cuts_found){
596 #if 0
597 new_cut->type = SUBTOUR_ELIM_SIDE;
598 new_cut->rhs = RHS(set_size, set_demand, capacity);
599 cuts_found += cg_send_cut(new_cut);
600 #endif
601 #ifdef DIRECTED_X_VARS
602 SEND_DIR_SUBTOUR_CONSTRAINT(set_size, set_demand);
603 #else
604 SEND_SUBTOUR_CONSTRAINT(set_size, set_demand);
605 #endif
606 memcpy(cutpt, coef, new_cut->size);
607 }
608
609 if ( cuts_found > max_num_cuts){
610 FREE(new_cut->coef);
611 return(cuts_found);
612 }
613 }
614 for (maxval = -1, pt = in_set+begin, dpt = cut_val+begin,
615 j = begin; j < end; pt++, dpt++, j++){
616 if (!(*pt) && *dpt > maxval){
617 maxval = cut_val[j];
618 max_vert = j;
619 }
620 }
621 if (maxval > 0){ /* add the vertex to the set */
622 in_set[max_vert]=1;
623 set_size+=1+ verts[compmembers[max_vert]].orig_node_list_size;
624 set_demand += demand[compmembers[max_vert]];
625 set_weight += verts[compmembers[max_vert]].weight;
626 cut_val[max_vert]=0;
627 for (e = verts[compmembers[max_vert]].first; e;
628 e = e->next_edge){
629 other_end = ref[e->other_end];
630 weight = e->data->weight;
631 set_weight += (in_set[other_end]) ? weight : -weight;
632 cut_val[other_end]+=(in_set[other_end]) ? 0 : weight;
633 }
634 }
635 else{ /* can't add anything to the set */
636 break;
637 }
638 }
639 }
640 }
641
642 FREE(new_cut->coef);
643 return(cuts_found);
644 }
645
646 /*===========================================================================*/
647
648 #if defined(ADD_FLOW_VARS) && defined(DIRECTED_X_VARS)
greedy_shrinking6_dicut(network * n,double capacity,double etol,cut_data * new_cut,int * compnodes,int * compmembers,int compnum,char * in_set,double * cut_val,int * ref,char * cut_list,int max_num_cuts,double * demand,int trial_num,double prob,int mult,int * num_cuts,int * alloc_cuts,cut_data *** cuts)649 int greedy_shrinking6_dicut(network *n, double capacity, double etol,
650 cut_data *new_cut, int *compnodes,
651 int *compmembers, int compnum,char *in_set,
652 double *cut_val, int *ref, char *cut_list,
653 int max_num_cuts, double *demand, int trial_num,
654 double prob, int mult, int *num_cuts,
655 int *alloc_cuts, cut_data ***cuts)
656 {
657 double set_demand, set_cut_val, tmp_cut_val;
658 vertex *verts = n->verts;
659 elist *e;
660 int cuts_found = 0, i, j;
661 char *pt, *cutpt;
662 double *dpt;
663 int vertnum = n->vertnum;
664
665 int min_vert = 0, set_size, begin = 1, cur_comp, end = 1, num_trials;
666 char *coef;
667 double min_val;
668 int max_size, numarcs, *arcs;
669 double denominator=pow(2.0,31.0)-1.0;
670 double r, q;
671
672 max_size = DSIZE + ISIZE + (vertnum >> DELETE_POWER) + 1
673 + 2*vertnum*vertnum*ISIZE;
674 new_cut->coef = (char *) calloc(max_size, sizeof(char));
675 coef = new_cut->coef + DSIZE + ISIZE;
676 arcs = (int *) (coef + (vertnum >> DELETE_POWER) + 1);
677
678 *in_set=0;
679
680 for(i = 1; i < vertnum; i++){
681 if (verts[compmembers[i]].deleted) compmembers[i] = 0;
682 ref[compmembers[i]] = i;
683 }
684 *ref = 0;
685 /* ref is a reference array for compmembers: gives a place
686 in which a vertex is listed in compmembers */
687
688 for (cur_comp = 1; cur_comp <= compnum; begin += compnodes[cur_comp],
689 cur_comp++){
690 /* for every component */
691 if (compnodes[cur_comp] <= 7) continue;
692
693 for (num_trials = 0; num_trials < trial_num * compnodes[cur_comp];
694 num_trials++){
695 end = begin + compnodes[cur_comp];
696 /*initialize the data structures */
697 memset(in_set + begin, 0, compnodes[cur_comp] * sizeof(char));
698 set_size = 0;
699 set_demand = 0;
700 set_cut_val = 0;
701 for (i = begin; i < end; i++ ){
702 if (compmembers[i] == 0) continue;
703 /*__BEGIN_EXPERIMENTAL_SECTION__*/
704 #if 0
705 r = CCutil_lprand(&rand_state)/CC_PRANDMAX;
706 #endif
707 /*___END_EXPERIMENTAL_SECTION___*/
708 r = RANDOM()/denominator;
709 q = (prob/compnodes[cur_comp]);
710 if (r < q){
711 in_set[i] = 1;
712 set_size += 1 + verts[compmembers[i]].orig_node_list_size;
713 set_demand += demand[compmembers[i]];
714 for (e = verts[compmembers[i]].first; e; e = e-> next_edge){
715 if (in_set[ref[e->other_end]]){
716 if (compmembers[i] < e->other_end){
717 set_cut_val -=
718 MIN(e->data->flow1,
719 MIN(set_demand, capacity)*e->data->weight1);
720 }else{
721 set_cut_val -=
722 MIN(e->data->flow2,
723 MIN(set_demand, capacity)*e->data->weight2);
724 }
725 }else{
726 if (compmembers[i] < e->other_end){
727 set_cut_val +=
728 MIN(e->data->flow2,
729 MIN(set_demand, capacity)*e->data->weight2);
730 }else{
731 set_cut_val +=
732 MIN(e->data->flow1,
733 MIN(set_demand, capacity)*e->data->weight1);
734 }
735 }
736 }
737 }
738 }
739 while(TRUE){
740 if (set_cut_val + etol < set_demand){
741 memset(coef, 0, ((vertnum >> DELETE_POWER) + 1)*sizeof(char));
742 numarcs = 0;
743 for (j = begin; j < end; j++){
744 if (in_set[j]){
745 (coef[(compmembers[j]) >> DELETE_POWER]) |=
746 (1 << ((compmembers[j]) & DELETE_AND));
747 for (e = verts[compmembers[j]].first; e; e=e->next_edge){
748 if (e->other_end < compmembers[j]){
749 if (!in_set[ref[e->other_end]] &&
750 set_demand*e->data->weight1 >
751 e->data->flow1+etol){
752 arcs[numarcs << 1] = e->other_end;
753 arcs[(numarcs << 1) + 1] = compmembers[j];
754 numarcs++;
755 }
756 }else{
757 if (!in_set[ref[e->other_end]] &&
758 set_demand*e->data->weight2 >
759 e->data->flow2+etol){
760 arcs[numarcs << 1] = e->other_end;
761 arcs[(numarcs << 1) + 1] = compmembers[j];
762 numarcs++;
763 }
764 }
765 }
766 }
767 }
768 ((double *)(new_cut->coef))[0] = set_demand;
769 ((int *)(new_cut->coef + DSIZE))[0] = numarcs;
770 new_cut->size = DSIZE + ISIZE + (vertnum >> DELETE_POWER)
771 + 1 + 2 * numarcs * ISIZE;
772 new_cut->type = MIXED_DICUT;
773 new_cut->rhs = set_demand;
774 new_cut->name = CUT__SEND_TO_CP;
775 cuts_found += cg_send_cut(new_cut, num_cuts, alloc_cuts, cuts);
776 if (cuts_found > max_num_cuts){
777 FREE(new_cut->coef);
778 return(cuts_found);
779 }
780 }
781 for (min_val = 0, pt = in_set + begin, j = begin; j < end;
782 pt++, j++){
783 if (!in_set[j]){
784 for (tmp_cut_val = 0, e = verts[compmembers[j]].first;
785 e; e = e->next_edge){
786 if (in_set[ref[e->other_end]]){
787 if (compmembers[j] < e->other_end){
788 tmp_cut_val -=
789 MIN(e->data->flow1,
790 MIN(set_demand, capacity)*e->data->weight1);
791 }else{
792 tmp_cut_val -=
793 MIN(e->data->flow2,
794 MIN(set_demand, capacity)*e->data->weight2);
795 }
796 }else{
797 if (compmembers[j] < e->other_end){
798 tmp_cut_val +=
799 MIN(e->data->flow2,
800 MIN(set_demand, capacity)*e->data->weight2);
801 }else{
802 tmp_cut_val +=
803 MIN(e->data->flow1,
804 MIN(set_demand, capacity)*e->data->weight1);
805 }
806 }
807 }
808 if (tmp_cut_val < min_val - etol){
809 min_vert = j;
810 min_val = tmp_cut_val;
811 }
812 }
813 }
814 if (min_val < 0){ /* add the vertex to the set */
815 in_set[min_vert] = 1;
816 set_size += 1 +
817 verts[compmembers[min_vert]].orig_node_list_size;
818 set_demand += demand[compmembers[min_vert]];
819 set_cut_val += min_val;
820 }
821 else{ /* can't add anything to the set */
822 break;
823 }
824 }
825 }
826 }
827 FREE(new_cut->coef);
828 return(cuts_found);
829 }
830 #endif
831
832 /*===========================================================================*/
833
greedy_shrinking1_one(network * n,double capacity,double etol,int max_num_cuts,cut_data * new_cut,char * in_set,double * cut_val,char * cut_list,int num_routes,double * demand,int mult,int * num_cuts,int * alloc_cuts,cut_data *** cuts)834 int greedy_shrinking1_one(network *n, double capacity, double etol,
835 int max_num_cuts, cut_data *new_cut,char *in_set,
836 double *cut_val, char *cut_list, int num_routes,
837 double *demand, int mult, int *num_cuts,
838 int *alloc_cuts, cut_data ***cuts)
839 {
840
841 double set_weight, set_cut_val, set_demand;
842 vertex *verts = n->verts;
843 elist *e;
844 int i, j, k, cuts_found = 0;
845 char *pt, *cutpt;
846 double *dpt;
847 int vertnum = n->vertnum;
848 int max_vert = 0;
849 int set_size;
850 /* int flag=0; */
851
852 double complement_demand, total_demand = verts[0].demand;
853 double complement_cut_val;
854 int complement_size;
855 char *coef;
856 double maxval;
857 int other_end;
858 double weight;
859 vertex *cur_nodept;
860
861 new_cut->size = (vertnum >> DELETE_POWER) + 1;
862 new_cut->coef = coef = (char *) (calloc(new_cut->size,sizeof(char)));
863 memset(cut_list, 0, new_cut->size * (max_num_cuts + 1));
864
865 for (i = 1; i < vertnum; i++ ){
866 if (verts[i].deleted) continue;/* for every node as a starting one */
867 /*initialize the data structures */
868 memset(in_set, 0, vertnum*sizeof(char));
869 memset(cut_val, 0,vertnum* sizeof(double));
870 in_set[i] = 1;
871 set_size = 1 + verts[i].orig_node_list_size;
872 set_cut_val = 0;
873 set_weight = verts[i].weight;
874 for (e= verts[i].first; e; e = e-> next_edge){
875 cut_val[e->other_end] = e->data->weight;
876 set_cut_val += e->data->weight;
877 }
878 set_demand = demand[i];
879
880 while(TRUE){
881 if (set_weight > RHS(set_size, set_demand, capacity) + etol &&
882 set_size > 2){
883 memset(coef, 0, new_cut->size*sizeof(char));
884 /* printf("%d :", i); */
885 /* printf("%d ", j); */
886 for (j = 1; j < vertnum; j++)
887 if (in_set[j]){
888 cur_nodept = verts + j;
889 if (cur_nodept->orig_node_list_size)
890 for (k = 0; k < cur_nodept->orig_node_list_size; k++)
891 (coef[(cur_nodept->orig_node_list)[k] >>
892 DELETE_POWER]) |=
893 (1 << ((cur_nodept->orig_node_list)[k] &
894 DELETE_AND));
895 (coef[j>> DELETE_POWER]) |= (1 << ( j & DELETE_AND));
896 }
897 /* printf("%f ", set_demand);
898 printf("%f \n",set_cut_val);*/
899 for (k = 0, cutpt = cut_list; k < cuts_found; k++,
900 cutpt += new_cut->size)
901 if (!memcmp(coef, cutpt, new_cut->size*sizeof(char)))
902 break; /* same cuts */
903 if ( k >= cuts_found){
904 #if 0
905 new_cut->type = SUBTOUR_ELIM_SIDE;
906 new_cut->rhs = RHS(set_size, set_demand, capacity);
907 cuts_found += cg_send_cut(new_cut);
908 #endif
909 #ifdef DIRECTED_X_VARS
910 SEND_DIR_SUBTOUR_CONSTRAINT(set_size, set_demand);
911 #else
912 SEND_SUBTOUR_CONSTRAINT(set_size, set_demand);
913 #endif
914 memcpy(cutpt, coef, new_cut->size);
915 }
916
917 if ( cuts_found > max_num_cuts){
918 FREE(new_cut->coef);
919 return(cuts_found);
920 }
921 }
922 /* check the complement */
923
924 complement_demand = total_demand - set_demand;
925 complement_cut_val = set_cut_val- 2*(*cut_val) + 2*num_routes;
926 complement_size = vertnum - 1 - set_size;
927 if (complement_cut_val< mult*(ceil(complement_demand/capacity))-etol
928 && complement_size > 2){
929 memset(coef, 0, new_cut->size*sizeof(char));
930 for (j = 1; j < vertnum; j++){
931 if (!(in_set[j]) && !(verts[j].deleted)){
932 cur_nodept = verts + j;
933 if (cur_nodept->orig_node_list_size)
934 for (k = 0; k < cur_nodept->orig_node_list_size; k++)
935 (coef[(cur_nodept->orig_node_list)[k] >>
936 DELETE_POWER]) |=
937 (1 << ((cur_nodept->orig_node_list)[k] &
938 DELETE_AND));
939 (coef[j>> DELETE_POWER]) |= (1 << ( j & DELETE_AND));
940 }
941 }
942 for (k=0, cutpt = cut_list; k < cuts_found; k++,
943 cutpt += new_cut->size)
944 if (!memcmp(coef, cutpt, new_cut->size*sizeof(char))) break;
945 if ( k >= cuts_found){
946 #if 0
947 new_cut->type = SUBTOUR_ELIM_SIDE;
948 new_cut->rhs = RHS(complement_size, complement_demand,
949 capacity);
950 cuts_found += cg_send_cut(new_cut);
951 #endif
952 #ifdef DIRECTED_X_VARS
953 SEND_DIR_SUBTOUR_CONSTRAINT(complement_size, complement_demand);
954 #else
955 SEND_SUBTOUR_CONSTRAINT(complement_size, complement_demand);
956 #endif
957 memcpy(cutpt, coef, new_cut->size);
958 }
959
960 if (cuts_found > max_num_cuts){
961 FREE(new_cut->coef);
962 return(cuts_found);
963 }
964 }
965
966 for (maxval = -1, pt = in_set, dpt = cut_val,pt++, dpt++,
967 j = 1; j < vertnum; pt++, dpt++, j++){
968 if (!(*pt) && *dpt > maxval){
969 maxval = cut_val[j];
970 max_vert = j;
971 }
972 }
973 if (maxval > 0){ /* add the vertex to the set */
974 in_set[max_vert] = 1;
975 set_size += 1 + verts[max_vert].orig_node_list_size ;
976 set_demand += demand[max_vert];
977 set_weight += verts[max_vert].weight;
978 cut_val[max_vert] = 0;
979 for (e = verts[max_vert].first; e; e = e-> next_edge){
980 other_end = e->other_end;
981 weight = e->data->weight;
982 set_weight += (in_set[other_end]) ? weight : -weight;
983 set_cut_val += (in_set[other_end]) ? (-weight) : weight;
984 cut_val[other_end] += weight;
985 }
986 }
987 else{ /* can't add anything to the set */
988 break;
989 }
990 }
991 }
992 FREE(new_cut->coef);
993 return(cuts_found);
994 }
995
996 /*===========================================================================*/
997
greedy_shrinking6_one(network * n,double capacity,double etol,cut_data * new_cut,char * in_set,double * cut_val,int num_routes,char * cut_list,int max_num_cuts,double * demand,int trial_num,double prob,int mult,int * num_cuts,int * alloc_cuts,cut_data *** cuts)998 int greedy_shrinking6_one(network *n, double capacity,
999 double etol, cut_data *new_cut,
1000 char *in_set, double *cut_val, int num_routes,
1001 char *cut_list, int max_num_cuts, double *demand,
1002 int trial_num, double prob, int mult, int *num_cuts,
1003 int *alloc_cuts, cut_data ***cuts)
1004 {
1005
1006 double set_weight, set_cut_val, set_demand;
1007 vertex *verts=n->verts;
1008 elist *e;
1009 int i, j, k, cuts_found = 0;
1010 char *pt, *cutpt;
1011 double *dpt;
1012 int vertnum = n->vertnum;
1013
1014 int max_vert = 0, set_size, begin = 1, end = 1, num_trials;
1015 char *coef;
1016 double maxval, r, q;
1017 double denominator = pow(2.0, 31.0) - 1.0;
1018
1019 int other_end;
1020 double weight;
1021
1022 double complement_demand, total_demand = verts[0].demand;
1023 double complement_cut_val;
1024 int complement_size;
1025 vertex *cur_nodept;
1026 /* int flag=0;*/
1027
1028 new_cut->size = (vertnum >> DELETE_POWER) + 1;
1029 new_cut->coef = coef = (char *) (calloc(new_cut->size,sizeof(char)));
1030 memset(cut_list, 0, new_cut->size * (max_num_cuts +1));
1031
1032 *in_set = 0;
1033
1034 for (num_trials = 0; num_trials < trial_num*vertnum ; num_trials++){
1035
1036 /*initialize the data structures */
1037 memset(in_set, 0, vertnum*sizeof(char));
1038 memset(cut_val, 0,vertnum* sizeof(double));
1039
1040 set_cut_val = 0;
1041 set_size = 0;
1042 set_demand = 0;
1043 set_weight = 0;
1044 for (i = 1 ; i < vertnum; i++ ){
1045 if (verts[i].deleted) continue;
1046 /*__BEGIN_EXPERIMENTAL_SECTION__*/
1047 #if 0
1048 r = CCutil_lprand(&rand_state)/CC_PRANDMAX;
1049 #endif
1050 /*___END_EXPERIMENTAL_SECTION___*/
1051 r = RANDOM()/denominator;
1052 q = (prob/vertnum);
1053 if (r < q){
1054 in_set[i] = 1;
1055 set_size += 1 + verts[i].orig_node_list_size;
1056 set_demand += demand[i];
1057 set_weight += verts[i].weight;
1058 for (e = verts[i].first; e; e = e-> next_edge){
1059 other_end = e->other_end;
1060 weight = e->data->weight;
1061 set_weight += (in_set[other_end]) ? weight : -weight;
1062 set_cut_val += (in_set[other_end]) ? (-weight) : weight;
1063 cut_val[other_end] += (in_set[other_end]) ? 0 : weight;
1064 }
1065 }
1066 }
1067 while(set_size){
1068 if (set_weight > RHS(set_size, set_demand, capacity) + etol &&
1069 set_size > 2){
1070 memset(coef, 0, new_cut->size*sizeof(char));
1071 for (j = 1; j < vertnum; j++ ){
1072 if (in_set[j]){
1073 cur_nodept = verts + j;
1074 if (cur_nodept->orig_node_list_size)
1075 for (k = 0; k < cur_nodept->orig_node_list_size; k++)
1076 (coef[(cur_nodept->orig_node_list)[k] >>
1077 DELETE_POWER]) |=
1078 (1 << ((cur_nodept->orig_node_list)[k] &
1079 DELETE_AND));
1080 (coef[j>> DELETE_POWER]) |= (1 << ( j & DELETE_AND));
1081 }
1082 }
1083 for (k = 0, cutpt = cut_list; k < cuts_found; k++,
1084 cutpt += new_cut->size)
1085 if (!memcmp(coef, cutpt, new_cut->size*sizeof(char))) break;
1086 if ( k >= cuts_found){
1087 #if 0
1088 new_cut->type = SUBTOUR_ELIM_SIDE;
1089 new_cut->rhs = RHS(set_size, set_demand,
1090 capacity);
1091 cuts_found += cg_send_cut(new_cut);
1092 #endif
1093 #ifdef DIRECTED_X_VARS
1094 SEND_DIR_SUBTOUR_CONSTRAINT(set_size, set_demand);
1095 #else
1096 SEND_SUBTOUR_CONSTRAINT(set_size, set_demand);
1097 #endif
1098 memcpy(cutpt, coef, new_cut->size);
1099 }
1100 if (cuts_found > max_num_cuts){
1101 FREE(new_cut->coef);
1102 return(cuts_found);
1103 }
1104 }
1105
1106 /* check the complement */
1107
1108 complement_demand = total_demand - set_demand;
1109 complement_cut_val = set_cut_val - 2*(*cut_val) + 2*num_routes;
1110 complement_size = vertnum - 1 - set_size;
1111 if (complement_cut_val< mult*(ceil(complement_demand/capacity))-
1112 etol && complement_size > 2){
1113 memset(coef, 0, new_cut->size*sizeof(char));
1114 for (j = 1; j < vertnum; j++){
1115 if (!(in_set[j])&& !(verts[j].deleted)){
1116 cur_nodept = verts + j;
1117 if (cur_nodept->orig_node_list_size)
1118 for (k = 0; k < cur_nodept->orig_node_list_size; k++)
1119 (coef[(cur_nodept->orig_node_list)[k] >>
1120 DELETE_POWER]) |=
1121 (1 << ((cur_nodept->orig_node_list)[k] &
1122 DELETE_AND));
1123 (coef[j>> DELETE_POWER]) |= (1 << ( j & DELETE_AND));
1124 }
1125 }
1126 for (k = 0, cutpt = cut_list; k < cuts_found; k++,
1127 cutpt += new_cut->size)
1128 if (!memcmp(coef, cutpt, new_cut->size*sizeof(char))) break;
1129 if ( k >= cuts_found){
1130 #if 0
1131 new_cut->type = SUBTOUR_ELIM_SIDE;
1132 new_cut->rhs = RHS(complement_size, complement_demand,
1133 capacity);
1134 cuts_found += cg_send_cut(new_cut);
1135 #endif
1136 #ifdef DIRECTED_X_VARS
1137 SEND_DIR_SUBTOUR_CONSTRAINT(complement_size, complement_demand);
1138 #else
1139 SEND_SUBTOUR_CONSTRAINT(complement_size, complement_demand);
1140 #endif
1141 memcpy(cutpt, coef, new_cut->size);
1142 }
1143
1144 if (cuts_found > max_num_cuts){
1145 FREE(new_cut->coef);
1146 return(cuts_found);
1147 }
1148 }
1149
1150 for (maxval = -1, pt = in_set + begin, dpt = cut_val + begin,
1151 j = begin; j < end; pt++, dpt++, j++){
1152 if (!(*pt) && *dpt > maxval){
1153 maxval = cut_val[j];
1154 max_vert = j;
1155 }
1156 }
1157 if (maxval > 0){ /* add the vertex to the set */
1158 in_set[max_vert] = 1;
1159 set_size += 1 + verts[max_vert].orig_node_list_size ;
1160 set_demand += demand[max_vert];
1161 set_weight += verts[max_vert].weight;
1162 cut_val[max_vert] = 0;
1163 for (e = verts[max_vert].first; e; e = e-> next_edge){
1164 other_end = e->other_end;
1165 weight = e->data->weight;
1166 set_weight += (in_set[other_end]) ? weight : -weight;
1167 set_cut_val += (in_set[other_end]) ? (-weight) : weight;
1168 cut_val[other_end] += weight;
1169 }
1170 }
1171 else{ /* can't add anything to the set */
1172 break;
1173 }
1174 }
1175 }
1176
1177 FREE(new_cut->coef);
1178 return(cuts_found);
1179 }
1180
1181 /*===========================================================================*/
1182
greedy_shrinking2_one(network * n,double capacity,double etol,cut_data * new_cut,char * in_set,double * cut_val,int num_routes,double * demand,int mult,int * num_cuts,int * alloc_cuts,cut_data *** cuts)1183 int greedy_shrinking2_one(network *n, double capacity,
1184 double etol, cut_data *new_cut,
1185 char *in_set, double *cut_val, int num_routes,
1186 double *demand, int mult, int *num_cuts,
1187 int *alloc_cuts, cut_data ***cuts)
1188 {
1189
1190 double set_cut_val, set_demand;
1191 vertex *verts = n->verts;
1192 elist *e, *cur_edge1, *cur_edge2;
1193 int j, k, cuts_found = 0;
1194 char *pt;
1195 double *dpt;
1196 int vertnum = n->vertnum;
1197
1198 int max_vert = 0, set_size, begin = 1, end = 1;
1199 char *coef;
1200 double maxval;
1201
1202 int other_end;
1203 double weight;
1204
1205 double complement_demand, total_demand = verts[0].demand;
1206 double complement_cut_val;
1207 int complement_size;
1208 vertex *cur_nodept;
1209
1210 new_cut->size = (vertnum >> DELETE_POWER) + 1;
1211 new_cut->coef =coef= (char *) (calloc(new_cut->size,sizeof(char)));
1212
1213 *in_set=0;
1214
1215 for (cur_edge1 = verts[0].first; cur_edge1;
1216 cur_edge1 = cur_edge1->next_edge){
1217 for (cur_edge2 = cur_edge1->next_edge; cur_edge2;
1218 cur_edge2 = cur_edge2->next_edge){
1219
1220 /*initialize the data structures */
1221 memset(in_set, 0, vertnum*sizeof(char));
1222 memset(cut_val, 0,vertnum* sizeof(double));
1223
1224 set_cut_val = 2;
1225 set_size = 2 + cur_edge1->other->orig_node_list_size +
1226 cur_edge2->other->orig_node_list_size;
1227 set_demand = demand[cur_edge1->other_end] +
1228 demand[cur_edge2->other_end];
1229 in_set[cur_edge1->other_end] = 1;
1230
1231 for (e = verts[cur_edge1->other_end].first; e; e = e-> next_edge){
1232 cut_val[e->other_end] += e->data->weight;
1233 }
1234
1235 in_set[cur_edge2->other_end] = 1;
1236 for (e = verts[cur_edge2->other_end].first; e; e = e-> next_edge){
1237 other_end = e->other_end;
1238 weight = e->data->weight;
1239 set_cut_val += (in_set[other_end]) ? (-weight) : weight;
1240 cut_val[other_end] += (in_set[other_end]) ? 0 : weight;
1241 }
1242 while(set_size){
1243 if (set_cut_val < mult*(ceil(set_demand/capacity)) - etol &&
1244 set_size > 2){
1245 memset(coef, 0, new_cut->size*sizeof(char));
1246 for (j = 1; j < vertnum; j++ ){
1247 if (in_set[j]){
1248 cur_nodept = verts + j;
1249 if (cur_nodept->orig_node_list_size)
1250 for (k = 0; k < cur_nodept->orig_node_list_size; k++)
1251 (coef[(cur_nodept->orig_node_list)[k] >>
1252 DELETE_POWER]) |=
1253 (1 << ((cur_nodept->orig_node_list)[k] &
1254 DELETE_AND));
1255 (coef[j >> DELETE_POWER]) |= (1 << (j & DELETE_AND));
1256 }
1257 }
1258 #if 0
1259 new_cut->type = SUBTOUR_ELIM_SIDE;
1260 new_cut->rhs = RHS(set_size, set_demand, capacity);
1261 cuts_found += cg_send_cut(new_cut);
1262 #endif
1263 #ifdef DIRECTED_X_VARS
1264 SEND_DIR_SUBTOUR_CONSTRAINT(set_size, set_demand);
1265 #else
1266 SEND_SUBTOUR_CONSTRAINT(set_size, set_demand);
1267 #endif
1268 }
1269
1270 /* check the complement */
1271
1272 complement_demand = total_demand - set_demand;
1273 complement_cut_val = set_cut_val - 2*(*cut_val) + 2*num_routes;
1274 complement_size = vertnum - 1 - set_size;
1275 if (complement_cut_val< mult*(ceil(complement_demand/capacity))-etol
1276 && complement_size > 2){
1277 memset(coef, 0, new_cut->size*sizeof(char));
1278 for (j = 1; j < vertnum; j++){
1279 if (!in_set[j]){
1280 cur_nodept=verts + j;
1281 if (cur_nodept->orig_node_list_size)
1282 for (k = 0; k < cur_nodept->orig_node_list_size; k++)
1283 (coef[(cur_nodept->orig_node_list)[k] >>
1284 DELETE_POWER]) |=
1285 (1 << ((cur_nodept->orig_node_list)[k] &
1286 DELETE_AND));
1287 (coef[j>> DELETE_POWER]) |= (1 << ( j & DELETE_AND));
1288 }
1289 }
1290 #if 0
1291 new_cut->type = SUBTOUR_ELIM_SIDE;
1292 new_cut->rhs = RHS(complement_size, complement_demand, capacity);
1293 cuts_found += cg_send_cut(new_cut);
1294 #endif
1295 #ifdef DIRECTED_X_VARS
1296 SEND_DIR_SUBTOUR_CONSTRAINT(complement_size, complement_demand);
1297 #else
1298 SEND_SUBTOUR_CONSTRAINT(complement_size, complement_demand);
1299 #endif
1300 SEND_SUBTOUR_CONSTRAINT(complement_size, complement_demand);
1301 }
1302
1303 for (maxval = -1, pt = in_set+begin, dpt = cut_val+begin,
1304 j = begin; j < end; pt++, dpt++, j++){
1305 if (!(*pt) && *dpt > maxval){
1306 maxval = cut_val[j];
1307 max_vert = j;
1308 }
1309 }
1310 if (maxval > 0){ /* add the vertex to the set */
1311 in_set[max_vert] = 1;
1312 set_size += 1 + verts[max_vert].orig_node_list_size ;
1313 set_demand += demand[max_vert];
1314 cut_val[max_vert] = 0;
1315 for (e = verts[max_vert].first; e; e = e-> next_edge){
1316 other_end = e->other_end;
1317 weight = e->data->weight;
1318 set_cut_val += (in_set[other_end]) ? (-weight) : weight;
1319 cut_val[other_end] += weight;
1320 }
1321 }
1322 else{ /* can't add anything to the set */
1323 break;
1324 }
1325 }
1326 }
1327 }
1328 FREE(new_cut->coef);
1329 return(cuts_found);
1330 }
1331