1 | /* $Id$ $Revision$ */ |
2 | /* vim:set shiftwidth=4 ts=8: */ |
3 | |
4 | /************************************************************************* |
5 | * Copyright (c) 2011 AT&T Intellectual Property |
6 | * All rights reserved. This program and the accompanying materials |
7 | * are made available under the terms of the Eclipse Public License v1.0 |
8 | * which accompanies this distribution, and is available at |
9 | * http://www.eclipse.org/legal/epl-v10.html |
10 | * |
11 | * Contributors: See CVS logs. Details at http://www.graphviz.org/ |
12 | *************************************************************************/ |
13 | |
14 | |
15 | /* |
16 | * Network Simplex Algorithm for Ranking Nodes of a DAG |
17 | */ |
18 | |
19 | #include "render.h" |
20 | #include <setjmp.h> |
21 | |
22 | static int init_graph(graph_t *); |
23 | static void dfs_cutval(node_t * v, edge_t * par); |
24 | static int dfs_range(node_t * v, edge_t * par, int low); |
25 | static int x_val(edge_t * e, node_t * v, int dir); |
26 | #ifdef DEBUG |
27 | static void check_cycles(graph_t * g); |
28 | #endif |
29 | |
30 | #define LENGTH(e) (ND_rank(aghead(e)) - ND_rank(agtail(e))) |
31 | #define SLACK(e) (LENGTH(e) - ED_minlen(e)) |
32 | #define SEQ(a,b,c) (((a) <= (b)) && ((b) <= (c))) |
33 | #define TREE_EDGE(e) (ED_tree_index(e) >= 0) |
34 | |
35 | static jmp_buf jbuf; |
36 | static graph_t *G; |
37 | static int N_nodes, N_edges; |
38 | static int Minrank, Maxrank; |
39 | static int S_i; /* search index for enter_edge */ |
40 | static int Search_size; |
41 | #define SEARCHSIZE 30 |
42 | static nlist_t Tree_node; |
43 | static elist Tree_edge; |
44 | |
45 | static void add_tree_edge(edge_t * e) |
46 | { |
47 | node_t *n; |
48 | //fprintf(stderr,"add tree edge %p %s ", (void*)e, agnameof(agtail(e))) ; fprintf(stderr,"%s\n", agnameof(aghead(e))) ; |
49 | if (TREE_EDGE(e)) { |
50 | agerr(AGERR, "add_tree_edge: missing tree edge\n" ); |
51 | longjmp (jbuf, 1); |
52 | } |
53 | ED_tree_index(e) = Tree_edge.size; |
54 | Tree_edge.list[Tree_edge.size++] = e; |
55 | if (ND_mark(agtail(e)) == FALSE) |
56 | Tree_node.list[Tree_node.size++] = agtail(e); |
57 | if (ND_mark(aghead(e)) == FALSE) |
58 | Tree_node.list[Tree_node.size++] = aghead(e); |
59 | n = agtail(e); |
60 | ND_mark(n) = TRUE; |
61 | ND_tree_out(n).list[ND_tree_out(n).size++] = e; |
62 | ND_tree_out(n).list[ND_tree_out(n).size] = NULL; |
63 | if (ND_out(n).list[ND_tree_out(n).size - 1] == 0) { |
64 | agerr(AGERR, "add_tree_edge: empty outedge list\n" ); |
65 | longjmp (jbuf, 1); |
66 | } |
67 | n = aghead(e); |
68 | ND_mark(n) = TRUE; |
69 | ND_tree_in(n).list[ND_tree_in(n).size++] = e; |
70 | ND_tree_in(n).list[ND_tree_in(n).size] = NULL; |
71 | if (ND_in(n).list[ND_tree_in(n).size - 1] == 0) { |
72 | agerr(AGERR, "add_tree_edge: empty inedge list\n" ); |
73 | longjmp (jbuf, 1); |
74 | } |
75 | } |
76 | |
77 | static void exchange_tree_edges(edge_t * e, edge_t * f) |
78 | { |
79 | int i, j; |
80 | node_t *n; |
81 | |
82 | ED_tree_index(f) = ED_tree_index(e); |
83 | Tree_edge.list[ED_tree_index(e)] = f; |
84 | ED_tree_index(e) = -1; |
85 | |
86 | n = agtail(e); |
87 | i = --(ND_tree_out(n).size); |
88 | for (j = 0; j <= i; j++) |
89 | if (ND_tree_out(n).list[j] == e) |
90 | break; |
91 | ND_tree_out(n).list[j] = ND_tree_out(n).list[i]; |
92 | ND_tree_out(n).list[i] = NULL; |
93 | n = aghead(e); |
94 | i = --(ND_tree_in(n).size); |
95 | for (j = 0; j <= i; j++) |
96 | if (ND_tree_in(n).list[j] == e) |
97 | break; |
98 | ND_tree_in(n).list[j] = ND_tree_in(n).list[i]; |
99 | ND_tree_in(n).list[i] = NULL; |
100 | |
101 | n = agtail(f); |
102 | ND_tree_out(n).list[ND_tree_out(n).size++] = f; |
103 | ND_tree_out(n).list[ND_tree_out(n).size] = NULL; |
104 | n = aghead(f); |
105 | ND_tree_in(n).list[ND_tree_in(n).size++] = f; |
106 | ND_tree_in(n).list[ND_tree_in(n).size] = NULL; |
107 | } |
108 | |
109 | static |
110 | void init_rank(void) |
111 | { |
112 | int i, ctr; |
113 | nodequeue *Q; |
114 | node_t *v; |
115 | edge_t *e; |
116 | |
117 | Q = new_queue(N_nodes); |
118 | ctr = 0; |
119 | |
120 | for (v = GD_nlist(G); v; v = ND_next(v)) { |
121 | if (ND_priority(v) == 0) |
122 | enqueue(Q, v); |
123 | } |
124 | |
125 | while ((v = dequeue(Q))) { |
126 | ND_rank(v) = 0; |
127 | ctr++; |
128 | for (i = 0; (e = ND_in(v).list[i]); i++) |
129 | ND_rank(v) = MAX(ND_rank(v), ND_rank(agtail(e)) + ED_minlen(e)); |
130 | for (i = 0; (e = ND_out(v).list[i]); i++) { |
131 | if (--(ND_priority(aghead(e))) <= 0) |
132 | enqueue(Q, aghead(e)); |
133 | } |
134 | } |
135 | if (ctr != N_nodes) { |
136 | agerr(AGERR, "trouble in init_rank\n" ); |
137 | for (v = GD_nlist(G); v; v = ND_next(v)) |
138 | if (ND_priority(v)) |
139 | agerr(AGPREV, "\t%s %d\n" , agnameof(v), ND_priority(v)); |
140 | } |
141 | free_queue(Q); |
142 | } |
143 | |
144 | static edge_t *leave_edge(void) |
145 | { |
146 | edge_t *f, *rv = NULL; |
147 | int j, cnt = 0; |
148 | |
149 | j = S_i; |
150 | while (S_i < Tree_edge.size) { |
151 | if (ED_cutvalue(f = Tree_edge.list[S_i]) < 0) { |
152 | if (rv) { |
153 | if (ED_cutvalue(rv) > ED_cutvalue(f)) |
154 | rv = f; |
155 | } else |
156 | rv = Tree_edge.list[S_i]; |
157 | if (++cnt >= Search_size) |
158 | return rv; |
159 | } |
160 | S_i++; |
161 | } |
162 | if (j > 0) { |
163 | S_i = 0; |
164 | while (S_i < j) { |
165 | if (ED_cutvalue(f = Tree_edge.list[S_i]) < 0) { |
166 | if (rv) { |
167 | if (ED_cutvalue(rv) > ED_cutvalue(f)) |
168 | rv = f; |
169 | } else |
170 | rv = Tree_edge.list[S_i]; |
171 | if (++cnt >= Search_size) |
172 | return rv; |
173 | } |
174 | S_i++; |
175 | } |
176 | } |
177 | return rv; |
178 | } |
179 | |
180 | static edge_t *Enter; |
181 | static int Low, Lim, Slack; |
182 | |
183 | static void dfs_enter_outedge(node_t * v) |
184 | { |
185 | int i, slack; |
186 | edge_t *e; |
187 | |
188 | for (i = 0; (e = ND_out(v).list[i]); i++) { |
189 | if (TREE_EDGE(e) == FALSE) { |
190 | if (!SEQ(Low, ND_lim(aghead(e)), Lim)) { |
191 | slack = SLACK(e); |
192 | if ((slack < Slack) || (Enter == NULL)) { |
193 | Enter = e; |
194 | Slack = slack; |
195 | } |
196 | } |
197 | } else if (ND_lim(aghead(e)) < ND_lim(v)) |
198 | dfs_enter_outedge(aghead(e)); |
199 | } |
200 | for (i = 0; (e = ND_tree_in(v).list[i]) && (Slack > 0); i++) |
201 | if (ND_lim(agtail(e)) < ND_lim(v)) |
202 | dfs_enter_outedge(agtail(e)); |
203 | } |
204 | |
205 | static void dfs_enter_inedge(node_t * v) |
206 | { |
207 | int i, slack; |
208 | edge_t *e; |
209 | |
210 | for (i = 0; (e = ND_in(v).list[i]); i++) { |
211 | if (TREE_EDGE(e) == FALSE) { |
212 | if (!SEQ(Low, ND_lim(agtail(e)), Lim)) { |
213 | slack = SLACK(e); |
214 | if ((slack < Slack) || (Enter == NULL)) { |
215 | Enter = e; |
216 | Slack = slack; |
217 | } |
218 | } |
219 | } else if (ND_lim(agtail(e)) < ND_lim(v)) |
220 | dfs_enter_inedge(agtail(e)); |
221 | } |
222 | for (i = 0; (e = ND_tree_out(v).list[i]) && (Slack > 0); i++) |
223 | if (ND_lim(aghead(e)) < ND_lim(v)) |
224 | dfs_enter_inedge(aghead(e)); |
225 | } |
226 | |
227 | static edge_t *enter_edge(edge_t * e) |
228 | { |
229 | node_t *v; |
230 | int outsearch; |
231 | |
232 | /* v is the down node */ |
233 | if (ND_lim(agtail(e)) < ND_lim(aghead(e))) { |
234 | v = agtail(e); |
235 | outsearch = FALSE; |
236 | } else { |
237 | v = aghead(e); |
238 | outsearch = TRUE; |
239 | } |
240 | Enter = NULL; |
241 | Slack = INT_MAX; |
242 | Low = ND_low(v); |
243 | Lim = ND_lim(v); |
244 | if (outsearch) |
245 | dfs_enter_outedge(v); |
246 | else |
247 | dfs_enter_inedge(v); |
248 | return Enter; |
249 | } |
250 | |
251 | static void init_cutvalues(void) |
252 | { |
253 | dfs_range(GD_nlist(G), NULL, 1); |
254 | dfs_cutval(GD_nlist(G), NULL); |
255 | } |
256 | |
257 | /* functions for initial tight tree construction */ |
258 | // borrow field from network simplex - overwritten in init_cutvalues() forgive me |
259 | #define ND_subtree(n) (subtree_t*)ND_par(n) |
260 | #define ND_subtree_set(n,value) (ND_par(n) = (edge_t*)value) |
261 | |
262 | typedef struct subtree_s { |
263 | node_t *rep; /* some node in the tree */ |
264 | int size; /* total tight tree size */ |
265 | int heap_index; /* required to find non-min elts when merged */ |
266 | struct subtree_s *par; /* union find */ |
267 | } subtree_t; |
268 | |
269 | /* find initial tight subtrees */ |
270 | static int tight_subtree_search(Agnode_t *v, subtree_t *st) |
271 | { |
272 | Agedge_t *e; |
273 | int i; |
274 | int rv; |
275 | |
276 | rv = 1; |
277 | ND_subtree_set(v,st); |
278 | for (i = 0; (e = ND_in(v).list[i]); i++) { |
279 | if (TREE_EDGE(e)) continue; |
280 | if ((ND_subtree(agtail(e)) == 0) && (SLACK(e) == 0)) { |
281 | add_tree_edge(e); |
282 | rv += tight_subtree_search(agtail(e),st); |
283 | } |
284 | } |
285 | for (i = 0; (e = ND_out(v).list[i]); i++) { |
286 | if (TREE_EDGE(e)) continue; |
287 | if ((ND_subtree(aghead(e)) == 0) && (SLACK(e) == 0)) { |
288 | add_tree_edge(e); |
289 | rv += tight_subtree_search(aghead(e),st); |
290 | } |
291 | } |
292 | return rv; |
293 | } |
294 | |
295 | static subtree_t *find_tight_subtree(Agnode_t *v) |
296 | { |
297 | subtree_t *rv; |
298 | rv = NEW(subtree_t); |
299 | rv->rep = v; |
300 | rv->size = tight_subtree_search(v,rv); |
301 | rv->par = rv; |
302 | return rv; |
303 | } |
304 | |
305 | typedef struct STheap_s { |
306 | subtree_t **elt; |
307 | int size; |
308 | } STheap_t; |
309 | |
310 | static subtree_t *STsetFind(Agnode_t *n0) |
311 | { |
312 | subtree_t *s0 = ND_subtree(n0); |
313 | while (s0->par && (s0->par != s0)) { |
314 | if (s0->par->par) {s0->par = s0->par->par;} /* path compression for the code weary */ |
315 | s0 = s0->par; |
316 | } |
317 | return s0; |
318 | } |
319 | |
320 | static subtree_t *STsetUnion(subtree_t *s0, subtree_t *s1) |
321 | { |
322 | subtree_t *r0, *r1, *r; |
323 | |
324 | for (r0 = s0; r0->par && (r0->par != r0); r0 = r0->par); |
325 | for (r1 = s1; r1->par && (r1->par != r1); r1 = r1->par); |
326 | if (r0 == r1) return r0; /* safety code but shouldn't happen */ |
327 | assert((r0->heap_index > -1) || (r1->heap_index > -1)); |
328 | if (r1->heap_index == -1) r = r0; |
329 | else if (r0->heap_index == -1) r = r1; |
330 | else if (r1->size < r0->size) r = r0; |
331 | else r = r1; |
332 | |
333 | r0->par = r1->par = r; |
334 | r->size = r0->size + r1->size; |
335 | assert(r->heap_index >= 0); |
336 | return r; |
337 | } |
338 | |
339 | #define INCIDENT(e,treeset) ((STsetFind(agtail(e),treeset)) != STsetFind(aghead(e),treeset)) |
340 | |
341 | /* find tightest edge to another tree incident on the given tree */ |
342 | static Agedge_t *inter_tree_edge_search(Agnode_t *v, Agnode_t *from, Agedge_t *best) |
343 | { |
344 | int i; |
345 | Agedge_t *e; |
346 | subtree_t *ts = STsetFind(v); |
347 | if (best && SLACK(best) == 0) return best; |
348 | for (i = 0; (e = ND_out(v).list[i]); i++) { |
349 | if (TREE_EDGE(e)) { |
350 | if (aghead(e) == from) continue; // do not search back in tree |
351 | best = inter_tree_edge_search(aghead(e),v,best); // search forward in tree |
352 | } |
353 | else { |
354 | if (STsetFind(aghead(e)) != ts) { // encountered candidate edge |
355 | if ((best == 0) || (SLACK(e) < SLACK(best))) best = e; |
356 | } |
357 | /* else ignore non-tree edge between nodes in the same tree */ |
358 | } |
359 | } |
360 | /* the following code must mirror the above, but for in-edges */ |
361 | for (i = 0; (e = ND_in(v).list[i]); i++) { |
362 | if (TREE_EDGE(e)) { |
363 | if (agtail(e) == from) continue; |
364 | best = inter_tree_edge_search(agtail(e),v,best); |
365 | } |
366 | else { |
367 | if (STsetFind(agtail(e)) != ts) { |
368 | if ((best == 0) || (SLACK(e) < SLACK(best))) best = e; |
369 | } |
370 | } |
371 | } |
372 | return best; |
373 | } |
374 | |
375 | static Agedge_t *inter_tree_edge(subtree_t *tree) |
376 | { |
377 | Agedge_t *rv; |
378 | rv = inter_tree_edge_search(tree->rep, (Agnode_t *)0, (Agedge_t *)0); |
379 | return rv; |
380 | } |
381 | |
382 | static |
383 | int STheapsize(STheap_t *heap) { return heap->size; } |
384 | |
385 | static |
386 | void STheapify(STheap_t *heap, int i) |
387 | { |
388 | int left, right, smallest; |
389 | subtree_t **elt = heap->elt; |
390 | do { |
391 | left = 2*(i+1)-1; |
392 | right = 2*(i+1); |
393 | if ((left < heap->size) && (elt[left]->size < elt[i]->size)) smallest = left; |
394 | else smallest = i; |
395 | if ((right < heap->size) && (elt[right]->size < elt[smallest]->size)) smallest = right; |
396 | else smallest = i; |
397 | if (smallest != i) { |
398 | subtree_t *temp; |
399 | temp = elt[i]; |
400 | elt[i] = elt[smallest]; |
401 | elt[smallest] = temp; |
402 | elt[i]->heap_index = i; |
403 | elt[smallest]->heap_index = smallest; |
404 | i = smallest; |
405 | } |
406 | else break; |
407 | } while (i < heap->size); |
408 | } |
409 | |
410 | static |
411 | STheap_t *STbuildheap(subtree_t **elt, int size) |
412 | { |
413 | int i; |
414 | STheap_t *heap; |
415 | heap = NEW(STheap_t); |
416 | heap->elt = elt; |
417 | heap->size = size; |
418 | for (i = 0; i < heap->size; i++) heap->elt[i]->heap_index = i; |
419 | for (i = heap->size/2; i >= 0; i--) |
420 | STheapify(heap,i); |
421 | return heap; |
422 | } |
423 | |
424 | static |
425 | subtree_t *(STheap_t *heap) |
426 | { |
427 | subtree_t *rv; |
428 | rv = heap->elt[0]; |
429 | rv->heap_index = -1; |
430 | heap->elt[0] = heap->elt[heap->size - 1]; |
431 | heap->elt[0]->heap_index = 0; |
432 | heap->elt[heap->size -1] = rv; /* needed to free storage later */ |
433 | heap->size--; |
434 | STheapify(heap,0); |
435 | return rv; |
436 | } |
437 | |
438 | static |
439 | void tree_adjust(Agnode_t *v, Agnode_t *from, int delta) |
440 | { |
441 | int i; |
442 | Agedge_t *e; |
443 | Agnode_t *w; |
444 | ND_rank(v) = ND_rank(v) + delta; |
445 | for (i = 0; (e = ND_tree_in(v).list[i]); i++) { |
446 | w = agtail(e); |
447 | if (w != from) |
448 | tree_adjust(w, v, delta); |
449 | } |
450 | for (i = 0; (e = ND_tree_out(v).list[i]); i++) { |
451 | w = aghead(e); |
452 | if (w != from) |
453 | tree_adjust(w, v, delta); |
454 | } |
455 | } |
456 | |
457 | static |
458 | subtree_t *merge_trees(Agedge_t *e) /* entering tree edge */ |
459 | { |
460 | int delta; |
461 | subtree_t *t0, *t1, *rv; |
462 | |
463 | assert(!TREE_EDGE(e)); |
464 | |
465 | t0 = STsetFind(agtail(e)); |
466 | t1 = STsetFind(aghead(e)); |
467 | |
468 | //fprintf(stderr,"merge trees of %d %d of %d, delta %d\n",t0->size,t1->size,N_nodes,delta); |
469 | |
470 | if (t0->heap_index == -1) { // move t0 |
471 | delta = SLACK(e); |
472 | tree_adjust(t0->rep,(Agnode_t*)0,delta); |
473 | } |
474 | else { // move t1 |
475 | delta = -SLACK(e); |
476 | tree_adjust(t1->rep,0,delta); |
477 | } |
478 | add_tree_edge(e); |
479 | rv = STsetUnion(t0,t1); |
480 | |
481 | return rv; |
482 | } |
483 | |
484 | /* Construct initial tight tree. Graph must be connected, feasible. |
485 | * Adjust ND_rank(v) as needed. add_tree_edge() on tight tree edges. |
486 | * trees are basically lists of nodes stored in nodequeues. |
487 | * Return 1 if input graph is not connected; 0 on success. |
488 | */ |
489 | static |
490 | int feasible_tree(void) |
491 | { |
492 | Agnode_t *n; |
493 | Agedge_t *ee; |
494 | subtree_t **tree, *tree0, *tree1; |
495 | int i, subtree_count = 0; |
496 | STheap_t *heap; |
497 | int error = 0; |
498 | |
499 | /* initialization */ |
500 | for (n = GD_nlist(G); n; n = ND_next(n)) { |
501 | ND_subtree_set(n,0); |
502 | } |
503 | |
504 | tree = N_NEW(N_nodes,subtree_t*); |
505 | /* given init_rank, find all tight subtrees */ |
506 | for (n = GD_nlist(G); n; n = ND_next(n)) { |
507 | if (ND_subtree(n) == 0) { |
508 | tree[subtree_count] = find_tight_subtree(n); |
509 | subtree_count++; |
510 | } |
511 | } |
512 | |
513 | /* incrementally merge subtrees */ |
514 | heap = STbuildheap(tree,subtree_count); |
515 | while (STheapsize(heap) > 1) { |
516 | tree0 = STextractmin(heap); |
517 | if (!(ee = inter_tree_edge(tree0))) { |
518 | error = 1; |
519 | break; |
520 | } |
521 | tree1 = merge_trees(ee); |
522 | STheapify(heap,tree1->heap_index); |
523 | } |
524 | |
525 | free(heap); |
526 | for (i = 0; i < subtree_count; i++) free(tree[i]); |
527 | free(tree); |
528 | if (error) return 1; |
529 | assert(Tree_edge.size == N_nodes - 1); |
530 | init_cutvalues(); |
531 | return 0; |
532 | } |
533 | |
534 | /* utility functions for debugging */ |
535 | static subtree_t *nd_subtree(Agnode_t *n) {return ND_subtree(n);} |
536 | static int nd_priority(Agnode_t *n) {return ND_priority(n);} |
537 | static int nd_rank(Agnode_t *n) {return ND_rank(n);} |
538 | static int ed_minlen(Agedge_t *e) {return ED_minlen(e);} |
539 | |
540 | /* walk up from v to LCA(v,w), setting new cutvalues. */ |
541 | static Agnode_t *treeupdate(Agnode_t * v, Agnode_t * w, int cutvalue, int dir) |
542 | { |
543 | edge_t *e; |
544 | int d; |
545 | |
546 | while (!SEQ(ND_low(v), ND_lim(w), ND_lim(v))) { |
547 | e = ND_par(v); |
548 | if (v == agtail(e)) |
549 | d = dir; |
550 | else |
551 | d = NOT(dir); |
552 | if (d) |
553 | ED_cutvalue(e) += cutvalue; |
554 | else |
555 | ED_cutvalue(e) -= cutvalue; |
556 | if (ND_lim(agtail(e)) > ND_lim(aghead(e))) |
557 | v = agtail(e); |
558 | else |
559 | v = aghead(e); |
560 | } |
561 | return v; |
562 | } |
563 | |
564 | static void rerank(Agnode_t * v, int delta) |
565 | { |
566 | int i; |
567 | edge_t *e; |
568 | |
569 | ND_rank(v) -= delta; |
570 | for (i = 0; (e = ND_tree_out(v).list[i]); i++) |
571 | if (e != ND_par(v)) |
572 | rerank(aghead(e), delta); |
573 | for (i = 0; (e = ND_tree_in(v).list[i]); i++) |
574 | if (e != ND_par(v)) |
575 | rerank(agtail(e), delta); |
576 | } |
577 | |
578 | /* e is the tree edge that is leaving and f is the nontree edge that |
579 | * is entering. compute new cut values, ranks, and exchange e and f. |
580 | */ |
581 | static void |
582 | update(edge_t * e, edge_t * f) |
583 | { |
584 | int cutvalue, delta; |
585 | Agnode_t *lca; |
586 | |
587 | delta = SLACK(f); |
588 | /* "for (v = in nodes in tail side of e) do ND_rank(v) -= delta;" */ |
589 | if (delta > 0) { |
590 | int s; |
591 | s = ND_tree_in(agtail(e)).size + ND_tree_out(agtail(e)).size; |
592 | if (s == 1) |
593 | rerank(agtail(e), delta); |
594 | else { |
595 | s = ND_tree_in(aghead(e)).size + ND_tree_out(aghead(e)).size; |
596 | if (s == 1) |
597 | rerank(aghead(e), -delta); |
598 | else { |
599 | if (ND_lim(agtail(e)) < ND_lim(aghead(e))) |
600 | rerank(agtail(e), delta); |
601 | else |
602 | rerank(aghead(e), -delta); |
603 | } |
604 | } |
605 | } |
606 | |
607 | cutvalue = ED_cutvalue(e); |
608 | lca = treeupdate(agtail(f), aghead(f), cutvalue, 1); |
609 | if (treeupdate(aghead(f), agtail(f), cutvalue, 0) != lca) { |
610 | agerr(AGERR, "update: mismatched lca in treeupdates\n" ); |
611 | longjmp (jbuf, 1); |
612 | } |
613 | ED_cutvalue(f) = -cutvalue; |
614 | ED_cutvalue(e) = 0; |
615 | exchange_tree_edges(e, f); |
616 | dfs_range(lca, ND_par(lca), ND_low(lca)); |
617 | } |
618 | |
619 | static void scan_and_normalize(void) |
620 | { |
621 | node_t *n; |
622 | |
623 | Minrank = INT_MAX; |
624 | Maxrank = -INT_MAX; |
625 | for (n = GD_nlist(G); n; n = ND_next(n)) { |
626 | if (ND_node_type(n) == NORMAL) { |
627 | Minrank = MIN(Minrank, ND_rank(n)); |
628 | Maxrank = MAX(Maxrank, ND_rank(n)); |
629 | } |
630 | } |
631 | if (Minrank != 0) { |
632 | for (n = GD_nlist(G); n; n = ND_next(n)) |
633 | ND_rank(n) -= Minrank; |
634 | Maxrank -= Minrank; |
635 | Minrank = 0; |
636 | } |
637 | } |
638 | |
639 | static void |
640 | freeTreeList (graph_t* g) |
641 | { |
642 | node_t *n; |
643 | for (n = GD_nlist(G); n; n = ND_next(n)) { |
644 | free_list(ND_tree_in(n)); |
645 | free_list(ND_tree_out(n)); |
646 | ND_mark(n) = FALSE; |
647 | } |
648 | } |
649 | |
650 | static void LR_balance(void) |
651 | { |
652 | int i, delta; |
653 | edge_t *e, *f; |
654 | |
655 | for (i = 0; i < Tree_edge.size; i++) { |
656 | e = Tree_edge.list[i]; |
657 | if (ED_cutvalue(e) == 0) { |
658 | f = enter_edge(e); |
659 | if (f == NULL) |
660 | continue; |
661 | delta = SLACK(f); |
662 | if (delta <= 1) |
663 | continue; |
664 | if (ND_lim(agtail(e)) < ND_lim(aghead(e))) |
665 | rerank(agtail(e), delta / 2); |
666 | else |
667 | rerank(aghead(e), -delta / 2); |
668 | } |
669 | } |
670 | freeTreeList (G); |
671 | } |
672 | |
673 | static int decreasingrankcmpf(node_t **n0, node_t **n1) { |
674 | return ND_rank(*n1) - ND_rank(*n0); |
675 | } |
676 | |
677 | static int increasingrankcmpf(node_t **n0, node_t **n1) { |
678 | return ND_rank(*n0) - ND_rank(*n1); |
679 | } |
680 | |
681 | static void TB_balance(void) |
682 | { |
683 | node_t *n; |
684 | edge_t *e; |
685 | int i, ii, low, high, choice, *nrank; |
686 | int inweight, outweight; |
687 | int adj = 0; |
688 | char *s; |
689 | |
690 | scan_and_normalize(); |
691 | |
692 | /* find nodes that are not tight and move to less populated ranks */ |
693 | nrank = N_NEW(Maxrank + 1, int); |
694 | for (i = 0; i <= Maxrank; i++) |
695 | nrank[i] = 0; |
696 | if ( (s = agget(G,"TBbalance" )) ) { |
697 | if (streq(s,"min" )) adj = 1; |
698 | else if (streq(s,"max" )) adj = 2; |
699 | if (adj) for (n = GD_nlist(G); n; n = ND_next(n)) |
700 | if (ND_node_type(n) == NORMAL) |
701 | if (ND_out(n).size == 0) |
702 | ND_rank(n) = ((adj == 1)? Minrank : Maxrank); |
703 | } |
704 | for (ii = 0, n = GD_nlist(G); n; ii++, n = ND_next(n)) { |
705 | Tree_node.list[ii] = n; |
706 | } |
707 | Tree_node.size = ii; |
708 | qsort(Tree_node.list, Tree_node.size, sizeof(Tree_node.list[0]), |
709 | adj > 1? decreasingrankcmpf : increasingrankcmpf); |
710 | for (i = 0; i < Tree_node.size; i++) { |
711 | n = Tree_node.list[i]; |
712 | if (ND_node_type(n) == NORMAL) |
713 | nrank[ND_rank(n)]++; |
714 | } |
715 | for (ii = 0; ii < Tree_node.size; ii++) { |
716 | n = Tree_node.list[ii]; |
717 | if (ND_node_type(n) != NORMAL) |
718 | continue; |
719 | inweight = outweight = 0; |
720 | low = 0; |
721 | high = Maxrank; |
722 | for (i = 0; (e = ND_in(n).list[i]); i++) { |
723 | inweight += ED_weight(e); |
724 | low = MAX(low, ND_rank(agtail(e)) + ED_minlen(e)); |
725 | } |
726 | for (i = 0; (e = ND_out(n).list[i]); i++) { |
727 | outweight += ED_weight(e); |
728 | high = MIN(high, ND_rank(aghead(e)) - ED_minlen(e)); |
729 | } |
730 | if (low < 0) |
731 | low = 0; /* vnodes can have ranks < 0 */ |
732 | if (adj) { |
733 | if (inweight == outweight) |
734 | ND_rank(n) = (adj == 1? low : high); |
735 | } |
736 | else { |
737 | if (inweight == outweight) { |
738 | choice = low; |
739 | for (i = low + 1; i <= high; i++) |
740 | if (nrank[i] < nrank[choice]) |
741 | choice = i; |
742 | nrank[ND_rank(n)]--; |
743 | nrank[choice]++; |
744 | ND_rank(n) = choice; |
745 | } |
746 | } |
747 | free_list(ND_tree_in(n)); |
748 | free_list(ND_tree_out(n)); |
749 | ND_mark(n) = FALSE; |
750 | } |
751 | free(nrank); |
752 | } |
753 | |
754 | static int init_graph(graph_t * g) |
755 | { |
756 | int i, feasible; |
757 | node_t *n; |
758 | edge_t *e; |
759 | |
760 | G = g; |
761 | N_nodes = N_edges = S_i = 0; |
762 | for (n = GD_nlist(g); n; n = ND_next(n)) { |
763 | ND_mark(n) = FALSE; |
764 | N_nodes++; |
765 | for (i = 0; (e = ND_out(n).list[i]); i++) |
766 | N_edges++; |
767 | } |
768 | |
769 | Tree_node.list = ALLOC(N_nodes, Tree_node.list, node_t *); |
770 | Tree_node.size = 0; |
771 | Tree_edge.list = ALLOC(N_nodes, Tree_edge.list, edge_t *); |
772 | Tree_edge.size = 0; |
773 | |
774 | feasible = TRUE; |
775 | for (n = GD_nlist(g); n; n = ND_next(n)) { |
776 | ND_priority(n) = 0; |
777 | for (i = 0; (e = ND_in(n).list[i]); i++) { |
778 | ND_priority(n)++; |
779 | ED_cutvalue(e) = 0; |
780 | ED_tree_index(e) = -1; |
781 | if (feasible |
782 | && (ND_rank(aghead(e)) - ND_rank(agtail(e)) < ED_minlen(e))) |
783 | feasible = FALSE; |
784 | } |
785 | ND_tree_in(n).list = N_NEW(i + 1, edge_t *); |
786 | ND_tree_in(n).size = 0; |
787 | for (i = 0; (e = ND_out(n).list[i]); i++); |
788 | ND_tree_out(n).list = N_NEW(i + 1, edge_t *); |
789 | ND_tree_out(n).size = 0; |
790 | } |
791 | return feasible; |
792 | } |
793 | |
794 | /* graphSize: |
795 | * Compute no. of nodes and edges in the graph |
796 | */ |
797 | static void |
798 | graphSize (graph_t * g, int* nn, int* ne) |
799 | { |
800 | int i, nnodes, nedges; |
801 | node_t *n; |
802 | edge_t *e; |
803 | |
804 | nnodes = nedges = 0; |
805 | for (n = GD_nlist(g); n; n = ND_next(n)) { |
806 | nnodes++; |
807 | for (i = 0; (e = ND_out(n).list[i]); i++) { |
808 | nedges++; |
809 | } |
810 | } |
811 | *nn = nnodes; |
812 | *ne = nedges; |
813 | } |
814 | |
815 | /* rank: |
816 | * Apply network simplex to rank the nodes in a graph. |
817 | * Uses ED_minlen as the internode constraint: if a->b with minlen=ml, |
818 | * rank b - rank a >= ml. |
819 | * Assumes the graph has the following additional structure: |
820 | * A list of all nodes, starting at GD_nlist, and linked using ND_next. |
821 | * Out and in edges lists stored in ND_out and ND_in, even if the node |
822 | * doesn't have any out or in edges. |
823 | * The node rank values are stored in ND_rank. |
824 | * Returns 0 if successful; returns 1 if the graph was not connected; |
825 | * returns 2 if something seriously wrong; |
826 | */ |
827 | int rank2(graph_t * g, int balance, int maxiter, int search_size) |
828 | { |
829 | int iter = 0, feasible; |
830 | char *ns = "network simplex: " ; |
831 | edge_t *e, *f; |
832 | |
833 | #ifdef DEBUG |
834 | check_cycles(g); |
835 | #endif |
836 | if (Verbose) { |
837 | int nn, ne; |
838 | graphSize (g, &nn, &ne); |
839 | fprintf(stderr, "%s %d nodes %d edges maxiter=%d balance=%d\n" , ns, |
840 | nn, ne, maxiter, balance); |
841 | start_timer(); |
842 | } |
843 | feasible = init_graph(g); |
844 | if (!feasible) |
845 | init_rank(); |
846 | if (maxiter <= 0) { |
847 | freeTreeList (g); |
848 | return 0; |
849 | } |
850 | |
851 | if (search_size >= 0) |
852 | Search_size = search_size; |
853 | else |
854 | Search_size = SEARCHSIZE; |
855 | |
856 | if (setjmp (jbuf)) { |
857 | return 2; |
858 | } |
859 | |
860 | if (feasible_tree()) { |
861 | freeTreeList (g); |
862 | return 1; |
863 | } |
864 | while ((e = leave_edge())) { |
865 | f = enter_edge(e); |
866 | update(e, f); |
867 | iter++; |
868 | if (Verbose && (iter % 100 == 0)) { |
869 | if (iter % 1000 == 100) |
870 | fputs(ns, stderr); |
871 | fprintf(stderr, "%d " , iter); |
872 | if (iter % 1000 == 0) |
873 | fputc('\n', stderr); |
874 | } |
875 | if (iter >= maxiter) |
876 | break; |
877 | } |
878 | switch (balance) { |
879 | case 1: |
880 | TB_balance(); |
881 | break; |
882 | case 2: |
883 | LR_balance(); |
884 | break; |
885 | default: |
886 | scan_and_normalize(); |
887 | freeTreeList (G); |
888 | break; |
889 | } |
890 | if (Verbose) { |
891 | if (iter >= 100) |
892 | fputc('\n', stderr); |
893 | fprintf(stderr, "%s%d nodes %d edges %d iter %.2f sec\n" , |
894 | ns, N_nodes, N_edges, iter, elapsed_sec()); |
895 | } |
896 | return 0; |
897 | } |
898 | |
899 | int rank(graph_t * g, int balance, int maxiter) |
900 | { |
901 | char *s; |
902 | int search_size; |
903 | |
904 | if ((s = agget(g, "searchsize" ))) |
905 | search_size = atoi(s); |
906 | else |
907 | search_size = SEARCHSIZE; |
908 | |
909 | return rank2 (g, balance, maxiter, search_size); |
910 | } |
911 | |
912 | /* set cut value of f, assuming values of edges on one side were already set */ |
913 | static void x_cutval(edge_t * f) |
914 | { |
915 | node_t *v; |
916 | edge_t *e; |
917 | int i, sum, dir; |
918 | |
919 | /* set v to the node on the side of the edge already searched */ |
920 | if (ND_par(agtail(f)) == f) { |
921 | v = agtail(f); |
922 | dir = 1; |
923 | } else { |
924 | v = aghead(f); |
925 | dir = -1; |
926 | } |
927 | |
928 | sum = 0; |
929 | for (i = 0; (e = ND_out(v).list[i]); i++) |
930 | sum += x_val(e, v, dir); |
931 | for (i = 0; (e = ND_in(v).list[i]); i++) |
932 | sum += x_val(e, v, dir); |
933 | ED_cutvalue(f) = sum; |
934 | } |
935 | |
936 | static int x_val(edge_t * e, node_t * v, int dir) |
937 | { |
938 | node_t *other; |
939 | int d, rv, f; |
940 | |
941 | if (agtail(e) == v) |
942 | other = aghead(e); |
943 | else |
944 | other = agtail(e); |
945 | if (!(SEQ(ND_low(v), ND_lim(other), ND_lim(v)))) { |
946 | f = 1; |
947 | rv = ED_weight(e); |
948 | } else { |
949 | f = 0; |
950 | if (TREE_EDGE(e)) |
951 | rv = ED_cutvalue(e); |
952 | else |
953 | rv = 0; |
954 | rv -= ED_weight(e); |
955 | } |
956 | if (dir > 0) { |
957 | if (aghead(e) == v) |
958 | d = 1; |
959 | else |
960 | d = -1; |
961 | } else { |
962 | if (agtail(e) == v) |
963 | d = 1; |
964 | else |
965 | d = -1; |
966 | } |
967 | if (f) |
968 | d = -d; |
969 | if (d < 0) |
970 | rv = -rv; |
971 | return rv; |
972 | } |
973 | |
974 | static void dfs_cutval(node_t * v, edge_t * par) |
975 | { |
976 | int i; |
977 | edge_t *e; |
978 | |
979 | for (i = 0; (e = ND_tree_out(v).list[i]); i++) |
980 | if (e != par) |
981 | dfs_cutval(aghead(e), e); |
982 | for (i = 0; (e = ND_tree_in(v).list[i]); i++) |
983 | if (e != par) |
984 | dfs_cutval(agtail(e), e); |
985 | if (par) |
986 | x_cutval(par); |
987 | } |
988 | |
989 | static int dfs_range(node_t * v, edge_t * par, int low) |
990 | { |
991 | edge_t *e; |
992 | int i, lim; |
993 | |
994 | lim = low; |
995 | ND_par(v) = par; |
996 | ND_low(v) = low; |
997 | for (i = 0; (e = ND_tree_out(v).list[i]); i++) |
998 | if (e != par) |
999 | lim = dfs_range(aghead(e), e, lim); |
1000 | for (i = 0; (e = ND_tree_in(v).list[i]); i++) |
1001 | if (e != par) |
1002 | lim = dfs_range(agtail(e), e, lim); |
1003 | ND_lim(v) = lim; |
1004 | return lim + 1; |
1005 | } |
1006 | |
1007 | #ifdef DEBUG |
1008 | void tchk(void) |
1009 | { |
1010 | int i, n_cnt, e_cnt; |
1011 | node_t *n; |
1012 | edge_t *e; |
1013 | |
1014 | n_cnt = 0; |
1015 | e_cnt = 0; |
1016 | for (n = agfstnode(G); n; n = agnxtnode(G, n)) { |
1017 | n_cnt++; |
1018 | for (i = 0; (e = ND_tree_out(n).list[i]); i++) { |
1019 | e_cnt++; |
1020 | if (SLACK(e) > 0) |
1021 | fprintf(stderr, "not a tight tree %p" , e); |
1022 | } |
1023 | } |
1024 | if ((n_cnt != Tree_node.size) || (e_cnt != Tree_edge.size)) |
1025 | fprintf(stderr, "something missing\n" ); |
1026 | } |
1027 | |
1028 | void check_cutvalues(void) |
1029 | { |
1030 | node_t *v; |
1031 | edge_t *e; |
1032 | int i, save; |
1033 | |
1034 | for (v = agfstnode(G); v; v = agnxtnode(G, v)) { |
1035 | for (i = 0; (e = ND_tree_out(v).list[i]); i++) { |
1036 | save = ED_cutvalue(e); |
1037 | x_cutval(e); |
1038 | if (save != ED_cutvalue(e)) |
1039 | abort(); |
1040 | } |
1041 | } |
1042 | } |
1043 | |
1044 | int check_ranks(void) |
1045 | { |
1046 | int cost = 0; |
1047 | node_t *n; |
1048 | edge_t *e; |
1049 | |
1050 | for (n = agfstnode(G); n; n = agnxtnode(G, n)) { |
1051 | for (e = agfstout(G, n); e; e = agnxtout(G, e)) { |
1052 | cost += (ED_weight(e)) * abs(LENGTH(e)); |
1053 | if (ND_rank(aghead(e)) - ND_rank(agtail(e)) - ED_minlen(e) < 0) |
1054 | abort(); |
1055 | } |
1056 | } |
1057 | fprintf(stderr, "rank cost %d\n" , cost); |
1058 | return cost; |
1059 | } |
1060 | |
1061 | void checktree(void) |
1062 | { |
1063 | int i, n = 0, m = 0; |
1064 | node_t *v; |
1065 | edge_t *e; |
1066 | |
1067 | for (v = agfstnode(G); v; v = agnxtnode(G, v)) { |
1068 | for (i = 0; (e = ND_tree_out(v).list[i]); i++) |
1069 | n++; |
1070 | if (i != ND_tree_out(v).size) |
1071 | abort(); |
1072 | for (i = 0; (e = ND_tree_in(v).list[i]); i++) |
1073 | m++; |
1074 | if (i != ND_tree_in(v).size) |
1075 | abort(); |
1076 | } |
1077 | fprintf(stderr, "%d %d %d\n" , Tree_edge.size, n, m); |
1078 | } |
1079 | |
1080 | void check_fast_node(node_t * n) |
1081 | { |
1082 | node_t *nptr; |
1083 | nptr = GD_nlist(agraphof(n)); |
1084 | while (nptr && nptr != n) |
1085 | nptr = ND_next(nptr); |
1086 | assert(nptr != NULL); |
1087 | } |
1088 | |
1089 | static char* dump_node (node_t* n) |
1090 | { |
1091 | static char buf[50]; |
1092 | |
1093 | if (ND_node_type(n)) { |
1094 | sprintf(buf, "%p" , n); |
1095 | return buf; |
1096 | } |
1097 | else |
1098 | return agnameof(n); |
1099 | } |
1100 | |
1101 | static void dump_graph (graph_t* g) |
1102 | { |
1103 | int i; |
1104 | edge_t *e; |
1105 | node_t *n,*w; |
1106 | FILE* fp = fopen ("ns.gv" , "w" ); |
1107 | fprintf (fp, "digraph \"%s\" {\n" , agnameof(g)); |
1108 | for (n = GD_nlist(g); n; n = ND_next(n)) { |
1109 | fprintf (fp, " \"%s\"\n" , dump_node(n)); |
1110 | } |
1111 | for (n = GD_nlist(g); n; n = ND_next(n)) { |
1112 | for (i = 0; (e = ND_out(n).list[i]); i++) { |
1113 | fprintf (fp, " \"%s\"" , dump_node(n)); |
1114 | w = aghead(e); |
1115 | fprintf (fp, " -> \"%s\"\n" , dump_node(w)); |
1116 | } |
1117 | } |
1118 | |
1119 | fprintf (fp, "}\n" ); |
1120 | fclose (fp); |
1121 | } |
1122 | |
1123 | static node_t *checkdfs(graph_t* g, node_t * n) |
1124 | { |
1125 | edge_t *e; |
1126 | node_t *w,*x; |
1127 | int i; |
1128 | |
1129 | if (ND_mark(n)) |
1130 | return 0; |
1131 | ND_mark(n) = TRUE; |
1132 | ND_onstack(n) = TRUE; |
1133 | for (i = 0; (e = ND_out(n).list[i]); i++) { |
1134 | w = aghead(e); |
1135 | if (ND_onstack(w)) { |
1136 | dump_graph (g); |
1137 | fprintf(stderr, "cycle: last edge %lx %s(%lx) %s(%lx)\n" , |
1138 | (uint64_t)e, |
1139 | agnameof(n), (uint64_t)n, |
1140 | agnameof(w), (uint64_t)w); |
1141 | return w; |
1142 | } |
1143 | else { |
1144 | if (ND_mark(w) == FALSE) { |
1145 | x = checkdfs(g, w); |
1146 | if (x) { |
1147 | fprintf(stderr,"unwind %lx %s(%lx)\n" , |
1148 | (uint64_t)e, |
1149 | agnameof(n), (uint64_t)n); |
1150 | if (x != n) return x; |
1151 | fprintf(stderr,"unwound to root\n" ); |
1152 | fflush(stderr); |
1153 | abort(); |
1154 | return 0; |
1155 | } |
1156 | } |
1157 | } |
1158 | } |
1159 | ND_onstack(n) = FALSE; |
1160 | return 0; |
1161 | } |
1162 | |
1163 | void check_cycles(graph_t * g) |
1164 | { |
1165 | node_t *n; |
1166 | for (n = GD_nlist(g); n; n = ND_next(n)) |
1167 | ND_mark(n) = ND_onstack(n) = FALSE; |
1168 | for (n = GD_nlist(g); n; n = ND_next(n)) |
1169 | checkdfs(g, n); |
1170 | } |
1171 | #endif /* DEBUG */ |
1172 | |