1 #include <memory.h>
2 #include <stdlib.h>
3
4 #include "BB_constants.h"
5 #include "BB_macros.h"
6 #include "min_cut.h"
7 #include "vrp_macros.h"
8 #include "capforest.h"
9 #include "vrp_const.h"
10
11 /*---------------------------------------------------------------------------*\
12 | This routine finds small cuts in the network by applying the ideas of |
13 | Nagamochi and Ibaraki |
14 \*---------------------------------------------------------------------------*/
15
min_cut(cg_vrp_spec * vrp,network * n,double etol)16 int min_cut(cg_vrp_spec *vrp, network *n, double etol)
17 {
18 cg_user_params *user_par = &vrp->par;
19 float mincut = 2;
20 register elist *ne;
21 register vertex *nv = NULL;
22 register vertex *nu = NULL;
23 register vertex *nw = NULL;
24 int i, j;
25 float contr_above = 2;
26 float nodecut = 0;
27 int vertnum = n->vertnum;
28 int edgenum = n->edgenum;
29 char scannedmark = 0;
30 vertex **nen;
31 int itnum =0;
32 int cut_size = (vertnum >> DELETE_POWER) + 1;
33 int total_demand = 0, capacity = vrp->capacity, num_cuts = 0;
34 cut_data *cut;
35 char change_in_vertnum = TRUE, *coef;
36 float max_q = 0;
37 vertex *verts = n->verts;
38 int v0, v1, cur_edge1_end;
39 elist *cur_edge1, *cur_edge2, *cur_edge3, *prev_edge;
40
41 /*allocate memory for existing nodes list and the binary tree used by
42 capforest*/
43 nen = n->enodes = (vertex **) calloc (n->vertnum, sizeof(vertex *));
44 n->tnodes = (vertex **) calloc (n->vertnum, sizeof(vertex *));
45
46 /*set up the data structure for the vertices*/
47 for(i = 0, total_demand = 0; i < vertnum; i++){
48 verts[i].scanned = 0;
49 verts[i].orignodenum = i;
50 verts[i].orig_node_list_size = 1;
51 verts[i].orig_node_list = (char *) calloc (cut_size, sizeof(char));
52 verts[i].orig_node_list[i >> DELETE_POWER] |= (1 << (i & DELETE_AND));
53 total_demand += verts[i].demand;
54 }
55
56 /*shrink all the 1-edges in the current graph -- this operation preserves
57 all small cuts in the graph and reduces the size of the problem*/
58 num_cuts += shrink_one_edges(vrp, n, &vertnum, &edgenum, capacity, etol);
59
60 cut = (cut_data *) calloc (1, sizeof(cut_data));
61 cut->size = cut_size;
62 coef = (char *) calloc (cut_size, sizeof(char));
63
64 while (change_in_vertnum){/*if there hasn't been a change in the number of
65 vertices since the last iteration, then we are
66 done*/
67 itnum++;
68 change_in_vertnum = FALSE;
69
70 scannedmark = 1-scannedmark;
71 max_q = capforest(n, vertnum, scannedmark);/*construct the capacitated
72 forests
73 discussed by N-I*/
74 contr_above = MIN(contr_above, max_q);/*reset the contraction limit*/
75
76 /*--------------------------------------------------------------------*\
77 | run through the list of "super nodes" and contract any edge whose |
78 | weight is above the limit "contr_above" |
79 \*--------------------------------------------------------------------*/
80
81 for( i=1; i<vertnum; i++){
82 nv = nen[i]; /*get the next existing node*/
83 ne = nv->first; /*get the first edge in its adjacency list*/
84 v0 = nv->orignodenum;
85 if (!v0) continue; /*the depot should not be considered as part of any
86 cut*/
87
88 /*-------------------------------------------------------------------*\
89 | scan the adjacency list of this node and contract all eligible edges|
90 | in its adjacency list |
91 \*-------------------------------------------------------------------*/
92
93 while (ne != NULL){
94 if ((ne->data->q >= contr_above - etol) && (ne->data->v0) &&
95 (ne->data->v1)){
96 v1 = OTHER_END(ne, v0);
97 nw = verts+v1; /*nw is a pointer to the other end vertex of the
98 current node*/
99
100 /*-----------------------------------------------------------------*\
101 | merge super node v0 with super node v1. "orig_node_list" is a bit |
102 | array containing the list of the original indices of nodes that |
103 | have been contracted into this super node. the demand of the super|
104 | node is the sum of the demands of all the nodes that have been |
105 | merged into it. |
106 \*-----------------------------------------------------------------*/
107 for (j = 0; j <cut_size; j++)
108 nv->orig_node_list[j] |= nw->orig_node_list[j];
109 nv->orig_node_list_size+=nw->orig_node_list_size;
110 nv->demand += nw->demand;
111 nw->scanned = TRUE;
112
113 /*-----------------------------------------------------------------*\
114 | Now update the adjacency lists to reflect the contraction of edge |
115 | (v0, v1) To do this, we scan through the adjacency list of node v1|
116 | For each of these edges, we look at the other end vertex of the |
117 | edge Then there are two cases. if the other end vertex is not also|
118 | adjacent to vertex v0, then we simply change the edge data so that|
119 | it originates from v0 instead of v1. If on the other hand, the |
120 | vertex in question is also adjacent to v0, then there already |
121 | exists an edge from v0 to the vertex in question and in this case,|
122 | we have to sum the weights of the two edges that are now to |
123 | originate from v0 and merge them into one edge. |
124 \*-----------------------------------------------------------------*/
125 for (cur_edge1 = nw->first; cur_edge1;
126 cur_edge1 = cur_edge1->next_edge){
127 cur_edge1_end = OTHER_END(cur_edge1, v1);/*get the other end of
128 the edge*/
129
130 for (cur_edge2 = nv->first, prev_edge = NULL; cur_edge2;
131 prev_edge = cur_edge2, cur_edge2 = cur_edge2->next_edge){
132
133 /*---remove edge (v0, v1) from the adjacency list of v0--------*/
134 if ((OTHER_END(cur_edge2, v0)) == v1){
135 cur_edge2->data->weight = 0;
136 cur_edge2->data->deleted = TRUE;
137 if (prev_edge)
138 prev_edge->next_edge = cur_edge2->next_edge;
139 else
140 nv->first = cur_edge2->next_edge;
141 if (!cur_edge2->next_edge) nv->last = prev_edge;
142 nv->degree--;
143 if (!nv->first) break;
144 }
145
146 /*if cur_edge1_end is also adjacent to v0, then we must merge
147 two edges*/
148 else if ((OTHER_END(cur_edge2, v0)) == cur_edge1_end){
149 cur_edge2->data->weight += cur_edge1->data->weight;
150 cur_edge1->data->deleted = TRUE;
151 cur_edge1->data->weight = 0;
152 nu = verts+cur_edge1_end;
153
154 /*check whether we have found a violated cut*/
155 if ((user_par->do_extra_checking) && (cur_edge1_end) &&
156 (cur_edge2->data->weight > 2 -
157 BINS(nv->demand + nu->demand, capacity))){
158 cur_edge3 = nv->first;
159 nodecut = 0;
160 while (cur_edge3 != NULL){
161 if ((cur_edge3->data != cur_edge2->data) &&
162 (cur_edge3->data != ne->data))
163 nodecut += cur_edge3->data->weight;
164 cur_edge3 = cur_edge3->next_edge;
165 }
166 cur_edge3 = nu->first;
167 while (cur_edge3 != NULL){
168 if ((cur_edge3->data != cur_edge2->data) &&
169 (cur_edge3->data != cur_edge1->data))
170 nodecut += cur_edge3->data->weight;
171 cur_edge3 = cur_edge3->next_edge;
172 }
173 cur_edge3 = nw->first;
174 while (cur_edge3 != NULL){
175 if ((cur_edge3->data != cur_edge1->data) &&
176 (cur_edge3->data != ne->data))
177 nodecut += cur_edge3->data->weight;
178 cur_edge3 = cur_edge3->next_edge;
179 }
180 /*if the cut is violated, then send it back to the LP*/
181 if (nodecut <
182 2*(BINS(nv->demand + nu->demand, capacity))- etol){
183 cut->type =
184 (nv->orig_node_list_size+nu->orig_node_list_size <
185 n->vertnum/2 ? SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
186 cut->rhs =
187 (cut->type == SUBTOUR_ELIM_SIDE ?
188 RHS(nv->orig_node_list_size+nu->orig_node_list_size,
189 nv->demand+nu->demand, capacity):
190 2*(BINS(nv->demand + nu->demand, capacity)));
191 for (j = 0; j <cut_size; j++){
192 coef[j] = (nv->orig_node_list[j]|nu->orig_node_list[j]);
193 }
194 cut->coef = coef;
195 num_cuts += cg_send_cut(cut);
196 }
197 } /*if ((cur_edge1_end) && (cur_edge2->data->weight > 2 - ...*/
198
199 nu->degree--;
200
201 /*----update the adjacency list of cur_edge1_end------------*/
202 for (prev_edge = NULL, cur_edge3 = nu->first;;
203 prev_edge = cur_edge3, cur_edge3 = cur_edge3->next_edge){
204 if ((OTHER_END(cur_edge3, cur_edge1_end)) == v1){
205 if (prev_edge){
206 prev_edge->next_edge = cur_edge3->next_edge;
207 }
208 else{
209 nu->first = cur_edge3->next_edge;
210 }
211 if (!cur_edge3->next_edge){
212 nu->last = prev_edge;
213 }
214 break;
215 }
216 }
217 cur_edge1->data->weight = 0;
218 } /*--else if ((OTHER_END(cur_edge2, v0)) == cur_edge1_end)*/
219 } /*--for (cur_edge2 = nv->first, prev_edge = NULL; cur_edge2;...*/
220 if (!nv->first) break;
221 } /*--for (cur_edge1 = nw->first; cur_edge1;...*/
222
223 /*-----------------------------------------------------------------*\
224 | Now merge the remaining adjacency lists of v0 and v1 |
225 \*-----------------------------------------------------------------*/
226
227 for (cur_edge1 = nw->first; cur_edge1;
228 cur_edge1 = cur_edge1->next_edge){
229 if (cur_edge1){
230 if (cur_edge1->data->weight){
231 if (nv->last){
232 nv->last->next_edge = cur_edge1;
233 nv->last = cur_edge1;
234 nv->degree++;
235 if (cur_edge1->data->v0 == v1)
236 cur_edge1->data->v0 = v0;
237 else
238 cur_edge1->data->v1 = v0;
239 }
240 else{
241 nv->first = nv->last = cur_edge1;
242 if (cur_edge1->data->v0 == v1)
243 cur_edge1->data->v0 = v0;
244 else
245 cur_edge1->data->v1 = v0;
246 }
247 }
248 }
249 }
250 nv->last->next_edge = NULL;
251 nw->degree = 0;
252 nw->first = nw->last = NULL;
253
254 /*-----------------------------------------------------------------*\
255 | Now we update the existing node list. We want to delete v1 from |
256 | the list but we have to be careful when we do this. We want to |
257 | make sure that all the still existing nodes that have not been |
258 | scanned yet in this iteration of the algorithm remain in a |
259 | position in the array such that they will be scanned. In other |
260 | words, in one of the positions with index >= i. There are several |
261 | cases that must be onsidered depending on where the node to be |
262 | deleted is positioned in the existing nodes array |
263 \*-----------------------------------------------------------------*/
264
265 if (nw->enodenum == i-1){
266 nen[nw->enodenum] = nen[--vertnum];
267 nen[vertnum]->enodenum = nw->enodenum;
268 i--;
269 }
270 else if (nw->enodenum < i-1){
271 nen[nw->enodenum] = nen[i-1];
272 nen[i-1]->enodenum = nw->enodenum;
273 nen[i-1] = nen[--vertnum];
274 nen[vertnum]->enodenum = i-1;
275 i--;
276 }
277 else{
278 nen[nw->enodenum] = nen[--vertnum];
279 nen[vertnum]->enodenum = nw->enodenum;
280 }
281 change_in_vertnum = TRUE;
282 edgenum--;
283
284 /*now check to see if the cut induced by super node v0 is violated*/
285 cur_edge1 = nv->first;
286 nodecut = 0;
287 while (cur_edge1 != NULL){
288 nodecut += cur_edge1->data->weight;
289 cur_edge1 = cur_edge1->next_edge;
290 }
291 if (nodecut < 2*(BINS(nv->demand, capacity)) - etol){
292 cut->type = (nv->orig_node_list_size < n->vertnum/2 ?
293 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
294 cut->rhs = (cut->type == SUBTOUR_ELIM_SIDE ?
295 RHS(nv->orig_node_list_size, nv->demand, capacity):
296 2*(BINS(nv->demand, capacity)));
297 cut->coef = nv->orig_node_list;
298 num_cuts += cg_send_cut(cut);
299 }
300 if (nodecut < mincut-etol){
301 mincut = nodecut;
302 if (user_par->update_contr_above)
303 contr_above = mincut;
304 }
305 ne = ne->next_edge;
306 } /*---------if ((ne->data->q >= contr_above - etol) ...*/
307 else{
308 ne = ne->next_edge;
309 } /*---------if ((ne->data->q >= contr_above - etol) ...*/
310 } /*-----------while (ne != NULL) */
311 } /*-------------for( i=1; i<vertnum; i++)*/
312
313 /*-----------------------------------------------------------------------*\
314 | finally as an extra check, we can scan through all the nodes after each |
315 | iteration is completely finished and see if the cuts induced by any of |
316 | these nodes are violated but this probably isn't necessary and may |
317 | result in the same cut being imposed several times |
318 \*-----------------------------------------------------------------------*/
319
320 #if 0
321 if ((vertnum > 1) && (user_par->do_extra_checking)){
322 for ( i = 1; i<vertnum; i++ ){
323 nodecut = 0;
324 nv = nen[i];
325 ne = nv->first;
326 while (ne != NULL){
327 nodecut += ne->data->weight;
328 ne = ne->next_edge;
329 }
330 if (nodecut < 2*(BINS(nv->demand, capacity)) - etol){
331 cut->type = (nv->orig_node_list_size < n->vertnum/2 ?
332 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
333 cut->rhs = (cut->type == SUBTOUR_ELIM_SIDE ?
334 RHS(nv->orig_node_list_size, nv->demand, capacity):
335 2*(BINS(nv->demand, capacity)));
336 cut->coef = nv->orig_node_list;
337 num_cuts += cg_send_cut(cut);
338 }
339 if (nodecut < mincut-etol){
340 mincut = nodecut;
341 if (user_par->update_contr_above)
342 contr_above = mincut;
343 }
344 }
345 }
346 else
347 #endif
348 contr_above = mincut;
349 }
350
351 FREE(coef);
352 FREE(cut);
353
354 n->mincut = mincut;
355
356 for (i = 0; i < n->vertnum; i++)
357 if (n->verts[i].orig_node_list)
358 free((char *)n->verts[i].orig_node_list);
359 free((char *)n->tnodes);
360 free((char *)n->enodes);
361 return(num_cuts);
362 }
363
364 /*===========================================================================*/
365
366 /*---------------------------------------------------------------------------*\
367 | This function shrinks 1-edges in the graph that are not adjacent to the |
368 | depot |
369 \*---------------------------------------------------------------------------*/
370
shrink_one_edges(cg_vrp_spec * vrp,network * n,int * cur_verts,int * cur_edges,int capacity,double etol)371 int shrink_one_edges(cg_vrp_spec *vrp, network *n, int *cur_verts,
372 int *cur_edges, int capacity, double etol)
373 {
374 cg_user_params *user_par = &vrp->par;
375 register vertex *nw;
376 register vertex *nv = NULL;
377 register vertex *nu;
378 register elist *ne;
379 register edge *dat = NULL;
380 vertex **nen = n->enodes;
381 int cut_size = (n->vertnum >> DELETE_POWER) + 1;
382 vertex *verts = n->verts;
383 int num_cuts = 0, itnum = 0;
384 float nodecut;
385
386 elist *new, *ne2, *nep, *prev;
387 int vertnum = *cur_verts;
388 int edgenum = 0;
389 int i, j;
390 cut_data *cut;
391 char change_in_vertnum = TRUE, *coef;
392 int v0, v1, cur_edge1_end;
393 elist *cur_edge1, *cur_edge2, *cur_edge3, *np, *prev_edge;
394
395 /*----------------------------------------------------------------------*\
396 | First we contract all chains of 1-edges, called 1-paths. We simply look|
397 | at all nodes of degree 2 and keep following the 1-path in each |
398 | direction from that node until we reach a node that is not of degree 2 |
399 | Contracting chains is much easier to do than actually contracting all |
400 | the 1-edges and so we do it first |
401 \*----------------------------------------------------------------------*/
402
403 cut = (cut_data *) calloc (1, sizeof(cut_data));
404 cut->size = cut_size;
405 coef = (char *) calloc (cut_size, sizeof(char));
406
407 for (i = 1; i < vertnum - 1; i++){
408 /*-----------check whether we have a degree 2 node-----------*/
409 if ((verts[i].degree != 2) || (verts[i].scanned))
410 continue;
411
412 nv = verts+i;
413 ne2 = (ne = nv->first)->next_edge;
414 /*-----------------------------------------------------------------------*\
415 | follow the 1-path from vertex i until we hit either the depot or a node |
416 | that is not of degree 2 |
417 \*-----------------------------------------------------------------------*/
418 for (nw = ne->other; nw->degree==2 && nw->orignodenum != 0; nw=ne->other){
419 if (nw == nv){/*if we come back to the same node, that means we have a
420 subtour and we can simply impose the appropriate
421 constraint */
422 cut->coef = nv->orig_node_list;
423 cut->type = (nv->orig_node_list_size < n->vertnum/2 ?
424 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
425 cut->rhs = (cut->type == SUBTOUR_ELIM_SIDE ?
426 RHS(nv->orig_node_list_size, nv->demand, capacity):
427 2*(BINS(nv->demand, capacity)));
428 num_cuts += cg_send_cut(cut);
429 nv->scanned = TRUE;
430 nv->first = nv->last = NULL;
431 break;
432 }
433 dat = ne->data;
434 nep = ne;
435 ne = nw->first;
436 if (ne->data == dat)
437 ne = ne->next_edge;
438 if (!ne->other_end){
439 ne = nep;
440 break;
441 }
442
443 nw->scanned = TRUE;
444
445 /*As we go along, we contract all the nodes on the 1-path into node i*/
446 nv->demand += nw->demand;
447 for (j = 0; j <cut_size; j++)
448 nv->orig_node_list[j] |= nw->orig_node_list[j];
449 nv->orig_node_list_size+=nw->orig_node_list_size;
450 nw->first = nw->last = NULL;
451 if (nv->demand > capacity){
452 cut->coef = nv->orig_node_list;
453 cut->type = (nv->orig_node_list_size < n->vertnum/2 ?
454 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
455 cut->rhs = (cut->type == SUBTOUR_ELIM_SIDE ?
456 RHS(nv->orig_node_list_size, nv->demand, capacity):
457 2*(BINS(nv->demand, capacity)));
458 num_cuts += cg_send_cut(cut);
459 }
460 } /*--for (nw = ne->other; nw->degree==2 && nw->orignodenum != 0;
461 nw=ne->other)..*/
462
463 if (nv->scanned) continue; /*in this case we had a subtour and so we don't
464 need to go on*/
465 if (nw->orignodenum){ /*otherwise, nw is a pointer to one of the endpoints
466 of the 1-path. make new be the last edge on the
467 path that is adjacent to nw*/
468 dat = ne->data;
469 for (ne=nw->first; ne->data != dat; ne=ne->next_edge);
470 new = ne;
471 }
472 else{ /*this indicates that nw points to the depot which means that node i
473 is joined to the depot by a one edge. since we do not want to
474 contract one edges that are adjacent to the depot, we ignore this
475 edge and consider i to be the endpoint of the 1-path*/
476 if (!((ne = nv->first)->other_end))
477 ne=ne->next_edge;
478 new = ne;
479 }
480
481 /*-----------------------------------------------------------------------*\
482 | Now we find the other end point of the 1-path in a similar fashion |
483 \*-----------------------------------------------------------------------*/
484 ne = ne2;
485 for (nu = ne->other; nu->degree==2 && nu->orignodenum != 0;
486 nu = ne->other){
487 dat = ne->data;
488 nep = ne;
489 ne = nu->first;
490 if (ne->data == dat)
491 ne = ne->next_edge;
492 if (!ne->other_end){
493 ne = nep;
494 break;
495 }
496
497 nu->scanned = TRUE;
498
499 nv->demand += nu->demand;
500 for (j = 0; j <cut_size; j++)
501 nv->orig_node_list[j] |= nu->orig_node_list[j];
502 nv->orig_node_list_size+=nu->orig_node_list_size;
503 nu->first = nu->last = NULL;
504 if (nv->demand > capacity){
505 cut->coef = nv->orig_node_list;
506 cut->type = (nv->orig_node_list_size < n->vertnum/2 ?
507 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
508 cut->rhs = (cut->type == SUBTOUR_ELIM_SIDE ?
509 RHS(nv->orig_node_list_size, nv->demand, capacity):
510 2*(BINS(nv->demand, capacity)));
511 num_cuts += cg_send_cut(cut);
512 }
513 }
514 if(nu->orignodenum){
515 dat = ne->data;
516 for (ne=nu->first; ne->data != dat; ne=ne->next_edge);
517 }
518 else{
519 if (!((ne = nv->first)->other_end))
520 ne = ne->next_edge;
521 }
522
523 /*-------------------------------------------------------------------*\
524 | Now we update the adjacency lists appropriately. Remember that if |
525 | either nu or nw point to the depot, we must not consider these as |
526 | endpoints of the 1-path so we have several cases |
527 \*-------------------------------------------------------------------*/
528
529 if(nu->orignodenum && nw->orignodenum){
530 dat = ne->data = new->data;
531 ne->other = nw;
532 new->other = nu;
533 ne->data->v0 = nu->orignodenum;
534 ne->data->v1 = nw->orignodenum;
535 nv->scanned = TRUE;
536 /*------------- contract nv into nu ----------------------------------*/
537 nu->demand += nv->demand;
538 for (j = 0; j <cut_size; j++)
539 nu->orig_node_list[j] |= nv->orig_node_list[j];
540 nu->orig_node_list_size+=nv->orig_node_list_size;
541 nv->first = nv->last = NULL;
542 /*-------- if nu induces a violated constraint, then impose it -------*/
543 if (nu->demand > capacity){
544 cut->coef = nu->orig_node_list;
545 cut->type = (nu->orig_node_list_size < n->vertnum/2 ?
546 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
547 cut->rhs = (cut->type == SUBTOUR_ELIM_SIDE ?
548 RHS(nu->orig_node_list_size, nu->demand, capacity):
549 2*(BINS(nu->demand, capacity)));
550 num_cuts += cg_send_cut(cut);
551 }
552 }
553 else if (!nu->orignodenum){/*in this case, nu points to the depot*/
554 dat = ne->data = new->data;
555 ne->other = nw;
556 new->other = nv;
557 ne->data->v0 = nw->orignodenum;
558 ne->data->v1 = nv->orignodenum;
559 nu = nv;
560 }
561 else if (!nw->orignodenum){/*in this case, nw points to the depot*/
562 dat = ne->data = new->data;
563 ne->other = nv;
564 new->other = nu;
565 ne->data->v0 = nu->orignodenum;
566 ne->data->v1 = nv->orignodenum;
567 nw = nv;
568 }
569
570 /*--------------------------------------------------------------------*\
571 | finally we check to see if there is another edge connecting nu to nw.|
572 | If so, then there must be a violated constraint since this implies a |
573 | cut of less than two in the graph. Also, this requires a bit of care |
574 | in updating the adjacency lists since then we will get a duplicate |
575 | edge between nu and nw |
576 \*--------------------------------------------------------------------*/
577
578 for (prev = NULL, ne = nu->first; ne; prev = ne, ne = ne->next_edge){
579 if ((ne->other == nw) && (ne->data != dat)){
580 dat->weight += ne->data->weight;
581 ne->data->weight = 0;
582 ne->data->deleted = TRUE;
583 if (prev)
584 prev->next_edge = ne->next_edge;
585 else
586 nu->first = ne->next_edge;
587 if (!ne->next_edge)
588 nu->last = prev;
589 cut->type =
590 (nu->orig_node_list_size+nw->orig_node_list_size < n->vertnum/2 ?
591 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
592 cut->rhs = (cut->type == SUBTOUR_ELIM_SIDE ?
593 RHS(nu->orig_node_list_size+nw->orig_node_list_size,
594 nu->demand+nw->demand, capacity):
595 2*(BINS(nu->demand+nw->demand, capacity)));
596 for (j = 0; j <cut_size; j++){
597 coef[j] = nu->orig_node_list[j] | nw->orig_node_list[j];
598 }
599 cut->coef = coef;
600 num_cuts += cg_send_cut(cut);
601 for (prev = NULL, ne = nw->first; ne; prev = ne, ne = ne->next_edge){
602 if (ne->data->deleted){
603 if (prev)
604 prev->next_edge = ne->next_edge;
605 else
606 nw->first = ne->next_edge;
607 if(!ne->next_edge)
608 nw->last = prev;
609 break;
610 }
611 }
612 break;
613 }
614 }
615 }
616
617 /*----- here we construct the list of existing nodes ----------------*/
618 for (i=0, j = 0; i<vertnum; j++){
619 if (!verts[j].scanned){
620 nen[i] = verts+j;
621 verts[j].enodenum = i;
622 i++;
623 edgenum += verts[j].degree;
624 }else{
625 vertnum--;
626 }
627 }
628 edgenum /= 2;
629
630 if (!user_par->shrink_one_edges){
631 *cur_verts = vertnum;
632 *cur_edges = edgenum;
633 return(num_cuts);
634 }
635
636 /*-----------------------------------------------------------------------*\
637 | Here we can optionally shrink all the 1-edges in the graph. I am not |
638 | sure whether this is advantageous or not but it makes it easier to spot |
639 | violated cuts. These edges would get contracted anyway in the next phase|
640 | of this algorithm but it might make sense to contract them beforehand. |
641 | The rest of this function is uncommented because it is exactly the same |
642 | as the main loop in the min_cut routine except that we only consider |
643 | 1-edges for contraction instead of any edge above the threshold |
644 \*-----------------------------------------------------------------------*/
645
646 while (change_in_vertnum){
647 itnum++;
648 change_in_vertnum = FALSE;
649 for( i=1; i<vertnum; i++){
650 nv = nen[i];
651 ne = nv->first;
652 v0 = nv->orignodenum;
653 while (ne != NULL){
654 if ((ne->data->weight >= 1 - etol) && (ne->data->v0) &&
655 (ne->data->v1)){
656 v1 = OTHER_END(ne, v0);
657 nw = verts+v1;
658 for (j = 0; j <cut_size; j++)
659 nv->orig_node_list[j] |= nw->orig_node_list[j];
660 nv->orig_node_list_size+=nw->orig_node_list_size;
661 nv->demand += nw->demand;
662 nw->scanned = TRUE;
663 if (nv->demand > capacity){
664 cut->coef = nv->orig_node_list;
665 cut->type = (nv->orig_node_list_size < n->vertnum/2 ?
666 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
667 cut->rhs = (cut->type == SUBTOUR_ELIM_SIDE ?
668 RHS(nv->orig_node_list_size, nv->demand, capacity):
669 2*(BINS(nv->demand, capacity)));
670 num_cuts += cg_send_cut(cut);
671 }
672
673 for (cur_edge1 = nw->first; cur_edge1;
674 cur_edge1 = cur_edge1->next_edge){
675 if ((cur_edge1_end = OTHER_END(cur_edge1, v1)) == v0){
676 cur_edge1->data->weight = 0;
677 cur_edge1->data->deleted = TRUE;
678 }
679
680 for (cur_edge2 = nv->first, prev_edge = NULL; cur_edge2;
681 prev_edge = cur_edge2, cur_edge2 = cur_edge2->next_edge){
682 if ((OTHER_END(cur_edge2, v0)) == v1){
683 cur_edge2->data->weight = 0;
684 cur_edge2->data->deleted = TRUE;
685 if (prev_edge)
686 prev_edge->next_edge = cur_edge2->next_edge;
687 else
688 nv->first = cur_edge2->next_edge;
689 if (!cur_edge2->next_edge) verts[v0].last = prev_edge;
690 nv->degree--;
691 if (!nv->first) break;
692 }
693 else if ((OTHER_END(cur_edge2, v0)) == cur_edge1_end){
694 cur_edge2->data->weight += cur_edge1->data->weight;
695 cur_edge1->data->deleted = TRUE;
696 cur_edge1->data->weight = 0;
697 nu = verts + cur_edge1_end;
698 if ((cur_edge1_end) && (cur_edge2->data->weight > 2 -
699 BINS(nv->demand + nu->demand, capacity) + etol)){
700 cut->type =
701 (nv->orig_node_list_size+nu->orig_node_list_size <
702 n->vertnum/2 ? SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
703 cut->rhs =
704 (cut->type == SUBTOUR_ELIM_SIDE ?
705 RHS(nv->orig_node_list_size+nu->orig_node_list_size,
706 nv->demand+nu->demand, capacity):
707 2*(BINS(nv->demand+nu->demand, capacity)));
708 for (j = 0; j <cut_size; j++){
709 coef[j] = nv->orig_node_list[j] | nu->orig_node_list[j];
710 }
711 cut->coef = coef;
712 num_cuts += cg_send_cut(cut);
713 }
714 nu->degree--;
715 for (prev_edge = NULL, cur_edge3 = nu->first;;
716 prev_edge = cur_edge3, cur_edge3 = cur_edge3->next_edge){
717 if ((OTHER_END(cur_edge3, cur_edge1_end)) == v1){
718 if (prev_edge){
719 prev_edge->next_edge = cur_edge3->next_edge;
720 }
721 else{
722 nu->first = cur_edge3->next_edge;
723 }
724 if (!cur_edge3->next_edge){
725 nu->last = prev_edge;
726 }
727 break;
728 }
729 }
730 cur_edge1->data->weight = 0;
731 }
732 }
733 if (!nv->first) break;
734 }
735
736 for (cur_edge1 = nw->first; cur_edge1;
737 cur_edge1 = cur_edge1->next_edge){
738 if (cur_edge1){
739 if (cur_edge1->data->weight){
740 if (nv->last){
741 nv->last->next_edge = cur_edge1;
742 nv->last = cur_edge1;
743 nv->degree++;
744 if (cur_edge1->data->v0 == v1)
745 cur_edge1->data->v0 = v0;
746 else
747 cur_edge1->data->v1 = v0;
748 }
749 else{
750 nv->first = nv->last = cur_edge1;
751 if (cur_edge1->data->v0 == v1)
752 cur_edge1->data->v0 = v0;
753 else
754 cur_edge1->data->v1 = v0;
755 }
756 }
757 }
758 }
759 nv->last->next_edge = NULL;
760 nw->degree = 0;
761 nw->first = nw->last = NULL;
762 if (nw->enodenum == i-1){
763 nen[nw->enodenum] = nen[--vertnum];
764 nen[vertnum]->enodenum = nw->enodenum;
765 i--;
766 }
767 else if (nw->enodenum < i-1){
768 nen[nw->enodenum] = nen[i-1];
769 nen[i-1]->enodenum = nw->enodenum;
770 nen[i-1] = nen[--vertnum];
771 nen[vertnum]->enodenum = i-1;
772 i--;
773 }
774 else{
775 nen[nw->enodenum] = nen[--vertnum];
776 nen[vertnum]->enodenum = nw->enodenum;
777 }
778 change_in_vertnum = TRUE;
779 edgenum--;
780 np = ne;
781 ne = ne->next_edge;
782 }
783 else{
784 np = ne;
785 ne = ne->next_edge;
786 }
787 }
788 }
789 }
790
791 if ((vertnum > 1) && (user_par->do_extra_checking)){
792 for ( i = 1; i<vertnum; i++ ){
793 nodecut = 0;
794 nv = nen[i];
795 ne = nv->first;
796 while (ne != NULL){
797 nodecut += ne->data->weight;
798 ne = ne->next_edge;
799 }
800 if (nodecut < 2*(BINS(nv->demand, capacity)) - etol){
801 cut->type = (nv->orig_node_list_size < n->vertnum/2 ?
802 SUBTOUR_ELIM_SIDE:SUBTOUR_ELIM_ACROSS);
803 cut->rhs = (cut->type == SUBTOUR_ELIM_SIDE ?
804 RHS(nv->orig_node_list_size, nv->demand, capacity):
805 2*(BINS(nv->demand, capacity)));
806 cut->coef = nv->orig_node_list;
807 num_cuts += cg_send_cut(cut);
808 }
809 }
810 }
811
812 *cur_verts = vertnum;
813 *cur_edges = edgenum;
814
815 FREE(coef);
816 FREE(cut);
817
818 return(num_cuts);
819 }
820