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