1 /**************************************************************************/
2 /* */
3 /* EXPORTED FUNCTION: */
4 /* */
5 /* double flow (); */
6 /* */
7 /**************************************************************************/
8
9 #include "machdefs.h"
10 #include "Xsubtour.h"
11
12 #ifdef CC_PROTOTYPE_ANSI
13
14 static void
15 relabel (Xnode *),
16 setlabels (Xgraph *G, Xnode *, Xnode *),
17 backwards_bfs (Xnode *, int),
18 addtoq (Xnodeset *, Xnode *),
19 dfs2 (Xgraph *G, Xnode *start),
20 marknode (Xnode *n, int v);
21
22 static int
23 subtourcuts (Xgraph *G, Xcplane **list),
24 connectcuts (Xgraph *G, Xcplane **list),
25 dfs (Xgraph *G, Xnode *start);
26
27 static Xnode
28 *push (Xnode *, Xedge *),
29 *popfromq (Xnodeset *);
30
31 #else
32
33 static void
34 relabel (),
35 setlabels (),
36 backwards_bfs (),
37 addtoq (),
38 dfs2 (),
39 marknode ();
40
41 static int
42 subtourcuts (),
43 connectcuts (),
44 dfs ();
45
46 static Xnode
47 *push (),
48 *popfromq ();
49
50 #endif
51
52 #define INFINITY (1<<30)
53
54 #ifdef CC_PROTOTYPE_ANSI
Xflow(Xgraph * G,Xnode * s,Xnode * t,double bound)55 double Xflow (Xgraph *G, Xnode *s, Xnode *t, double bound)
56 #else
57 double Xflow (G, s, t, bound)
58 Xgraph *G;
59 Xnode *s, *t;
60 double bound;
61 #endif
62 {
63 Xnodeset q;
64 Xnode *n, *n1;
65 Xedge *e;
66 int savelabel;
67 Xedgeptr *ep;
68 int count, round;
69
70 q.head = q.tail = (Xnodeptr *) NULL;
71 for (n = G->pseudonodelist->next; n; n = n->next) {
72 n->excess = 0.0;
73 n->active = 0;
74 n->current = n->cadj.head;
75 n->flowlabel = 0;
76 for (ep = n->cadj.head; ep; ep = ep->next)
77 ep->this->flow = 0.0;
78 }
79 setlabels (G, s, t);
80 t->active = 1; /* a lie, which keeps t off the active queue */
81 for (ep = s->cadj.head; ep; ep = ep->next) {
82 e = ep->this;
83 if (!e->stay)
84 continue;
85 if (e->cends[0] == s) {
86 e->flow = e->x;
87 if (e->x > 0.0) {
88 e->cends[1]->excess += e->x;
89 addtoq (&q, e->cends[1]);
90 }
91 } else {
92 e->flow = -e->x;
93 if (e->x > 0.0) {
94 e->cends[0]->excess += e->x;
95 addtoq (&q, e->cends[0]);
96 }
97 }
98 }
99 count = 0;
100 round = G->npseudonodes / 2;
101 while (q.head && t->excess < bound) {
102 if (count == round) {
103 setlabels (G, s, t);
104 count = 0;
105 } else
106 count++;
107 n = popfromq (&q);
108 n->active = 0;
109 savelabel = n->flowlabel;
110 do {
111 ep = n->current;
112 if ((n1 = push (n, ep->this)) != (Xnode *) NULL)
113 addtoq (&q, n1);
114 else {
115 n->current = ep->next;
116 if (!n->current) {
117 n->current = n->cadj.head;
118 relabel (n);
119 }
120 }
121 } while (n->excess > 0.0 && n->flowlabel == savelabel);
122 if (n->excess > 0.0 && n->flowlabel < G->npseudonodes)
123 addtoq (&q, n);
124 }
125
126 while (q.head) {
127 popfromq (&q);
128 }
129
130 return t->excess;
131 }
132
133 #ifdef CC_PROTOTYPE_ANSI
relabel(Xnode * n)134 static void relabel (Xnode *n)
135 #else
136 static void relabel (n)
137 Xnode *n;
138 #endif
139 {
140 int m = INFINITY;
141 Xedgeptr *ep;
142 Xedge *e;
143 int t;
144
145 for (ep = n->cadj.head; ep; ep = ep->next) {
146 e = ep->this;
147 if (n == e->cends[0]) {
148 if (e->x - e->flow > 0.0 && (t = e->cends[1]->flowlabel) < m)
149 m = t;
150 } else {
151 if (e->x + e->flow > 0.0 && (t = e->cends[0]->flowlabel) < m)
152 m = t;
153 }
154 }
155 n->flowlabel = m + 1;
156 }
157
158 #ifdef CC_PROTOTYPE_ANSI
push(Xnode * n,Xedge * e)159 static Xnode *push (Xnode *n, Xedge *e)
160 #else
161 static Xnode *push (n, e)
162 Xnode *n;
163 Xedge *e;
164 #endif
165 {
166 Xnode *n1;
167 double rf;
168
169 if (e->cends[0] == n) {
170 n1 = e->cends[1];
171 rf = e->x - e->flow;
172 if (n->flowlabel == n1->flowlabel + 1 && rf > 0.0) {
173 if (n->excess < rf)
174 rf = n->excess;
175 n->excess -= rf;
176 e->flow += rf;
177 n1->excess += rf;
178 return n1;
179 } else
180 return 0;
181 } else {
182 n1 = e->cends[0];
183 rf = e->x + e->flow;
184 if (n->flowlabel == n1->flowlabel + 1 && rf > 0.0) {
185 if (n->excess < rf)
186 rf = n->excess;
187 n->excess -= rf;
188 e->flow -= rf;
189 n1->excess += rf;
190 return n1;
191 } else
192 return 0;
193 }
194 }
195
196 #ifdef CC_PROTOTYPE_ANSI
setlabels(Xgraph * G,Xnode * s,Xnode * t)197 static void setlabels (Xgraph *G, Xnode *s, Xnode *t)
198 #else
199 static void setlabels (G, s, t)
200 Xgraph *G;
201 Xnode *s, *t;
202 #endif
203 {
204 G->magicnum++;
205 t->flowlabel = 0;
206 backwards_bfs (t, G->magicnum);
207 s->flowlabel = G->npseudonodes;
208
209 if (s->magiclabel == G->magicnum)
210 return;
211 else
212 backwards_bfs (s, G->magicnum);
213 }
214
215 #ifdef CC_PROTOTYPE_ANSI
backwards_bfs(Xnode * s,int K)216 static void backwards_bfs (Xnode *s, int K)
217 #else
218 static void backwards_bfs (s, K)
219 Xnode *s;
220 int K;
221 #endif
222 {
223 Xnode *this, *next, *tail;
224 Xedge *e;
225 Xedgeptr *ep;
226 int dist;
227
228 s->magiclabel = K;
229 next = s;
230 s->tnext = (Xnode *) NULL;
231 do {
232 for (this = next, next = (Xnode *) NULL; this; this = this->tnext) {
233 dist = this->flowlabel + 1;
234 for (ep = this->cadj.head; ep; ep = ep->next) {
235 e = ep->this;
236 if (this == e->cends[0]) {
237 tail = e->cends[1];
238 if (tail->magiclabel != K &&
239 e->x + e->flow > 0.0) {
240 tail->flowlabel = dist;
241 tail->tnext = next;
242 next = tail;
243 tail->magiclabel = K;
244 }
245 } else {
246 tail = e->cends[0];
247 if (tail->magiclabel != K &&
248 e->x - e->flow > 0.0) {
249 tail->flowlabel = dist;
250 tail->tnext = next;
251 next = tail;
252 tail->magiclabel = K;
253 }
254 }
255 }
256 }
257 } while (next);
258 }
259
260 #ifdef CC_PROTOTYPE_ANSI
addtoq(Xnodeset * q,Xnode * n)261 static void addtoq (Xnodeset *q, Xnode *n)
262 #else
263 static void addtoq (q, n)
264 Xnodeset *q;
265 Xnode *n;
266 #endif
267 {
268 Xnodeptr *newn;
269
270 if (!n->active) {
271 newn = Xnodeptralloc ();
272 newn->next = (Xnodeptr *) NULL;
273 newn->this = n;
274 if (q->tail) {
275 q->tail->next = newn;
276 } else {
277 q->head = newn;
278 }
279 q->tail = newn;
280 n->active = 1;
281 }
282 }
283
284 #ifdef CC_PROTOTYPE_ANSI
popfromq(Xnodeset * q)285 static Xnode *popfromq (Xnodeset *q)
286 #else
287 static Xnode *popfromq (q)
288 Xnodeset *q;
289 #endif
290 {
291 Xnodeptr *newn = q->head;
292 Xnode *n = newn->this;
293
294 q->head = newn->next;
295 if (!q->head)
296 q->tail = (Xnodeptr *) NULL;
297 Xnodeptrfree (newn);
298 return n;
299 }
300
301 #ifdef CC_PROTOTYPE_ANSI
subtourcuts(Xgraph * G,Xcplane ** list)302 static int subtourcuts (Xgraph *G, Xcplane **list)
303 #else
304 static int subtourcuts (G, list)
305 Xgraph *G;
306 Xcplane **list;
307 #endif
308 {
309 Xnode *n, *nnext;
310 int returnval = 0;
311
312 for (n = G->pseudonodelist->next->next; n; n = nnext) {
313 /* printf ("o"); fflush (stdout); */
314
315 nnext = n->next;
316 if (Xflow (G, G->pseudonodelist->next, n, XCUTTWO) < XCUTTWO) {
317 G->magicnum++;
318 dfs2 (G, n);
319 returnval += Xloadcplane_cut (G, list, G->magicnum);
320 /* printf ("i"); fflush (stdout); */
321 Xsimpleshrink (G, G->pseudonodelist->next, n);
322 }
323 }
324 /* printf ("\n"); */
325 return returnval;
326 }
327
328 #ifdef CC_PROTOTYPE_ANSI
dfs2(Xgraph * G,Xnode * start)329 static void dfs2 (Xgraph *G, Xnode *start)
330 #else
331 static void dfs2 (G, start)
332 Xgraph *G;
333 Xnode *start;
334 #endif
335 {
336 Xedgeptr *epp;
337 Xedge *ep;
338 Xnodeptr *next, *queue = (Xnodeptr *) NULL;
339 Xnode *n;
340
341 marknode (start, G->magicnum);
342 Xadd_nodeptr (&queue, start);
343
344 while (queue) {
345 n = queue->this;
346 next = queue->next;
347 Xnodeptrfree (queue);
348 queue = next;
349 for (epp = n->cadj.head; epp; epp = epp->next) {
350 ep = epp->this;
351 if (ep->cends[0] == n) {
352 if (ep->x + ep->flow > 0.0 &&
353 ep->cends[1]->magiclabel != G->magicnum) {
354 marknode (ep->cends[1], G->magicnum);
355 Xadd_nodeptr (&queue, ep->cends[1]);
356 }
357 } else {
358 if (ep->x - ep->flow > 0.0 &&
359 ep->cends[0]->magiclabel != G->magicnum) {
360 marknode (ep->cends[0], G->magicnum);
361 Xadd_nodeptr (&queue, ep->cends[0]);
362 }
363 }
364 }
365 }
366 }
367
368 #ifdef CC_PROTOTYPE_ANSI
marknode(Xnode * n,int v)369 static void marknode (Xnode *n, int v)
370 #else
371 static void marknode (n, v)
372 Xnode *n;
373 int v;
374 #endif
375 {
376 Xnodeptr *np;
377
378 n->magiclabel = v;
379 for (np = n->base.head; np; np = np->next) {
380 np->this->magiclabel = v;
381 }
382 }
383
384 #ifdef CC_PROTOTYPE_ANSI
Xexactcutcheck(Xgraph * G,Xcplane ** list,double * x)385 int Xexactcutcheck (Xgraph *G, Xcplane **list, double *x)
386 #else
387 int Xexactcutcheck (G, list, x)
388 Xgraph *G;
389 Xcplane **list;
390 double *x;
391 #endif
392 {
393 int i, hit;
394
395 Xloadx (G, x);
396 Xbuildpseudonodelist (G, 0);
397 i = Xshrinkprocess (G, list);
398 /*
399 printf ("%d heavy-edge cut(s), np = %d\n", i, G->npseudonodes);
400 fflush (stdout);
401 */
402 if (i >= XCUTNUM) {
403 Xdestroypseudonodelist (G);
404 return i;
405 }
406 hit = i;
407 Xrebuildcadj (G);
408
409 i = subtourcuts (G, list);
410 /* printf ("%d flow cut(s)\n", i); fflush (stdout); */
411 Xdestroypseudonodelist (G);
412
413 return hit + i;
414 }
415
416 #ifdef CC_PROTOTYPE_ANSI
Xmincut(Xgraph * G,Xnode * s,Xnode * t,double bound,double * value,int * label)417 int Xmincut (Xgraph *G, Xnode *s, Xnode *t, double bound, double *value,
418 int *label)
419 #else
420 int Xmincut (G, s, t, bound, value, label)
421 Xgraph *G;
422 Xnode *s, *t;
423 double bound, *value;
424 int *label;
425 #endif
426 {
427 /* Uses the x field as capacities. If min cut from s to t has */
428 /* value < bound it returns 1 and marks nodes in cut with */
429 /* magiclabel equal to label. value get the capacity of cut. */
430
431 *value = Xflow (G, s, t, bound);
432 if (*value >= bound)
433 return 0;
434 else {
435 G->magicnum++;
436 dfs2 (G, t);
437 *label = G->magicnum;
438 return 1;
439 }
440 }
441
442 #ifdef CC_PROTOTYPE_ANSI
Xrunconnectcuts(Xgraph * G,Xcplane ** list,double * x)443 int Xrunconnectcuts (Xgraph *G, Xcplane **list, double *x)
444 #else
445 int Xrunconnectcuts (G, list, x)
446 Xgraph *G;
447 Xcplane **list;
448 double *x;
449 #endif
450 {
451 int i;
452
453 Xloadx (G, x);
454 Xbuildpseudonodelist (G, 0);
455 i = connectcuts (G, list);
456 Xdestroypseudonodelist (G);
457
458 fflush (stdout);
459 return i;
460 }
461
462 #ifdef CC_PROTOTYPE_ANSI
connectcuts(Xgraph * G,Xcplane ** list)463 static int connectcuts (Xgraph *G, Xcplane **list)
464 #else
465 static int connectcuts (G, list)
466 Xgraph *G;
467 Xcplane **list;
468 #endif
469 {
470 int val, ccount = 0, count, track = 1;
471 Xnode *n;
472
473 val = ++(G->magicnum);
474 if ((count = dfs (G, G->pseudonodelist->next)) < G->npseudonodes) {
475 ccount += Xloadcplane_cut (G, list, G->magicnum);
476 n = G->pseudonodelist->next;
477 do {
478 track++;
479 n = n->next;
480 while (n->magiclabel >= val)
481 n = n->next;
482 G->magicnum++;
483 count += dfs (G, n);
484 } while (count < G->npseudonodes);
485
486 if (track > 2) {
487 for (val++; val <= G->magicnum; val++) {
488 ccount += Xloadcplane_cut (G, list, val);
489 }
490 return ccount;
491 } else {
492 return ccount;
493 }
494 }
495 return 0;
496 }
497
498 #ifdef CC_PROTOTYPE_ANSI
dfs(Xgraph * G,Xnode * start)499 static int dfs (Xgraph *G, Xnode *start)
500 #else
501 static int dfs (G, start)
502 Xgraph *G;
503 Xnode *start;
504 #endif
505 {
506 int i = 0;
507 Xedgeptr *epp;
508 Xedge *ep;
509 Xnode *n, *n1;
510 Xnodeptr *next, *queue = (Xnodeptr *) NULL;
511
512 marknode (start, G->magicnum);
513 Xadd_nodeptr (&queue, start);
514
515 while (queue) {
516 i++;
517 n = queue->this;
518 next = queue->next;
519 Xnodeptrfree (queue);
520 queue = next;
521 for (epp = n->cadj.head; epp; epp = epp->next) {
522 ep = epp->this;
523 if (ep->stay) {
524 n1 = ep->cends[0];
525 if (n1 == n)
526 n1 = ep->cends[1];
527 if (n1->magiclabel != G->magicnum) {
528 marknode (n1, G->magicnum);
529 Xadd_nodeptr (&queue, n1);
530 }
531 }
532 }
533 }
534 return i;
535 }
536
537