![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dag_pos.c CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [Research Unix] / researchv10no / cmd / dag|
Return to dag_pos.c CVS log | Up to [Research Unix] / researchv10no / cmd / dag |
1.1 root 1: #include "draw_dag.h"
2:
3: #include <sys/types.h>
4: #include <sys/times.h>
5: #define TIC 60
6:
7:
8: /*
9: Assign horizontal positions to nodes
10: */
11:
12: // Edge factors depend on their adjacent nodes.
13: // The weight system chosen below encourages straight long edges and paths.
14: // The node naming scheme is:
15: // s: a real node with <=1 in-edge and <=1 out-edge
16: // v: a virtual node
17: // e: anything else.
18: static const int C_ee = 1,
19: C_es = 1,
20: C_ev = 1,
21: C_vs = 2,
22: C_ss = 2,
23: C_vv = 8;
24:
25: static int *Selfwidth; // width of nodes with self edges
26:
27: static void levelposition(int),
28: bestposition(),
29: goodposition(),
30: smoothedges();
31: static int rlength(int*, int*);
32: static inline int widthof(int);
33:
34:
35: #ifdef DEBUG
36: static void checkpos(int *nodepos)
37: {
38: for(int l = 0; l <= Maxlevel; ++l)
39: for(int *r = Rank[l]+1; *r != -1; ++r)
40: if((nodepos[r[-1]]+widthof(r[-1])/2+Nodesep) >
41: (nodepos[r[0]]-widthof(r[0])/2))
42: panic("Bad position assignment");
43: }
44: #endif
45:
46:
47: // ranksepar = 0 for full flexibility
48: // 1 for forcing all ranks be separate by a max amount
49: // computed to avoid edge/node intersections
50: // >1 for fixing at user's requested separation
51:
52: void dag_positions(int ranksepar, int opt_level)
53: {
54: struct tms begtm, endtm;
55: if(Verbose)
56: {
57: #ifndef TIMING
58: fprintf(stderr,"Assign positions\n");
59: #endif
60: times(&begtm);
61: }
62:
63: // additional widths for nodes with self edges
64: if(N_self_edges > 0)
65: {
66: Selfwidth = new int[N_nodes];
67: for(int i = 0; i < N_real_nodes; ++i)
68: {
69: if(Rank[Level[i]][Invert[i]+1] == -1)
70: continue;
71:
72: for(edge_t *e = Noedges[i]; e; e = e->next)
73: if(e->node == i)
74: break;
75: if(e)
76: Selfwidth[i] = e->width + Nodesep;
77: }
78: }
79:
80: // compute node positions
81: Nodepos = new int[N_nodes];
82: if(opt_level >= 1)
83: bestposition();
84: else goodposition();
85:
86: // readjust node positions if there are self edges
87: if(N_self_edges > 0)
88: {
89: for(int i = 0; i < N_real_nodes; ++i)
90: if(Selfwidth[i] > 0)
91: Nodepos[i] -= Selfwidth[i]/2;
92: delete Selfwidth;
93: }
94:
95: // compute layer position
96: Levelpos = new int[Maxlevel+1];
97: levelposition(ranksepar);
98:
99: // smooth out little turns in long edges
100: Deleted = new int[N_nodes];
101: smoothedges();
102:
103: // readjust positions to start from 0
104: int minx = Maxint;
105: for(int i = 0; i <= Maxlevel; ++i)
106: if(Nodepos[Rank[i][0]] < minx)
107: minx = Nodepos[Rank[i][0]];
108: if(minx != 0)
109: for(i = 0; i < N_nodes; ++i)
110: Nodepos[i] -= minx;
111:
112: if(Verbose)
113: {
114: times(&endtm);
115: int total = (endtm.tms_utime-begtm.tms_utime) +
116: (endtm.tms_stime-begtm.tms_stime);
117: int percent = (total - (total/TIC)*TIC)*100/TIC;
118: #ifdef TIMING
119: fprintf(stderr,"%d.%02d\t",total/TIC, percent);
120: #else
121: fprintf(stderr,"Time in coordinates assignment: %d.%02ds\n",
122: total/TIC, percent);
123: #endif
124: }
125: }
126:
127:
128:
129: // for an edge top->bot and the horizontal assignment of nodes,
130: // compute the separation between respective layers that would
131: // guarantee that the edge does not intersect any drawings of other
132: // nodes.
133:
134: static int vsepar(int top, int bot, int *toprank, int *botrank, int separ)
135: {
136: int xtop = Nodepos[top];
137: int xbot = Nodepos[bot];
138: int minx = min(xtop,xbot);
139: int maxx = max(xtop,xbot);
140: int htop = top < N_real_nodes ? Height[top]/2 : Levelheight[Level[top]]/2;
141: int hbot = bot < N_real_nodes ? Height[bot]/2 : Levelheight[Level[bot]]/2;
142:
143: if(maxx-minx < Nodesep)
144: return separ;
145:
146: for(; *toprank != -1; ++toprank)
147: if(Nodepos[*toprank] > minx && *toprank < N_real_nodes)
148: break;
149: for(; *toprank != -1; ++toprank)
150: {
151: if(Nodepos[*toprank] >= maxx)
152: break;
153: if(*toprank >= N_real_nodes)
154: continue;
155:
156: // the side of the box that may intersect with the edge
157: int x = Nodepos[*toprank] + (Width[*toprank]/2)*(xtop>xbot ? 1:-1);
158: if(x == xtop)
159: continue;
160:
161: // the height of that side
162: int h = Height[*toprank]+Nodesep;
163:
164: // the separation required
165: int s = (int)((double)((htop*(xbot-x)+h*(xtop-xbot)))/(xtop-x));
166: if((s += hbot) >= separ)
167: separ = s;
168: }
169: for(; *botrank != -1; ++botrank)
170: if(Nodepos[*botrank] > minx && *botrank < N_real_nodes)
171: break;
172: for(; *botrank != -1; ++botrank)
173: {
174: if(Nodepos[*botrank] >= maxx)
175: break;
176: if(*botrank >= N_real_nodes)
177: continue;
178: int x = Nodepos[*botrank] + (Width[*botrank]/2)*(xbot>xtop ? 1:-1);
179: if(x == xbot)
180: continue;
181: int h = Height[*botrank]+Nodesep;
182: int s = (int)((double)((hbot*(x-xtop)+h*(xtop-xbot)))/(x-xbot));
183: if((s += htop) > separ)
184: separ = s;
185: }
186:
187: return separ;
188: }
189:
190:
191:
192: static void levelposition(int ranksepar)
193: {
194: int maxsep = 0;
195: Levelpos[0] = Levelheight[0]/2;
196: for(int i = 1; i <= Maxlevel; ++i)
197: {
198: int h_min = Levelsep + (Levelheight[i-1]+Levelheight[i])/2;
199: int h_max = h_min+Levelsep;
200:
201: // reduce edge/node intersections by raising level separations
202: if(ranksepar < 2)
203: {
204: for(int *vp = Rank[i-1]; *vp != -1; ++vp)
205: for(edge_t *e = Outedges[*vp]; e; e = e->next)
206: h_min = vsepar(*vp,e->node,Rank[i-1],Rank[i],h_min);
207: if(ranksepar == 1)
208: maxsep = max(maxsep,min(h_max,h_min));
209: }
210:
211: if(ranksepar != 1)
212: Levelpos[i] = Levelpos[i-1] + min(h_max,h_min);
213: }
214:
215: if(ranksepar == 1)
216: for(i = 1; i <= Maxlevel; ++i)
217: Levelpos[i] = Levelpos[i-1]+maxsep;
218: }
219:
220: // see if a node is part of a path
221: static int pathdegree(int n)
222: {
223: int deg = 0;
224: if(Inedges[n])
225: {
226: deg += 1;
227: if(Inedges[n]->next)
228: deg += 3;
229: }
230: if(Outedges[n])
231: {
232: deg += 1;
233: if(Outedges[n]->next)
234: deg += 3;
235: }
236: return deg;
237: }
238:
239: // multiplicative factor for different type of edges
240: static int edgefactor(int v, int w)
241: {
242: #define V_BIT 001
243: #define S_BIT 002
244: #define C_EE 000 /* 0000 */
245: #define C_VE 004 /* 0100 */
246: #define C_SE 010 /* 1000 */
247: #define C_EV 001 /* 0001 */
248: #define C_VV 005 /* 0101 */
249: #define C_SV 011 /* 1001 */
250: #define C_ES 002 /* 0010 */
251: #define C_VS 006 /* 0110 */
252: #define C_SS 012 /* 1010 */
253:
254: int v_type = 0;
255: if(v >= N_real_nodes)
256: v_type |= V_BIT;
257: else if(pathdegree(v) <= 2)
258: v_type |= S_BIT;
259:
260: int w_type = 0;
261: if(w >= N_real_nodes)
262: w_type |= V_BIT;
263: else if(pathdegree(w) <= 2)
264: w_type |= S_BIT;
265:
266: switch((v_type<<2)|w_type)
267: {
268: default:
269: case C_EE:
270: return C_ee;
271: case C_VE:
272: return C_ev;
273: case C_EV:
274: return C_ev;
275: case C_SE:
276: return C_es;
277: case C_ES:
278: return C_es;
279: case C_VS:
280: return C_vs;
281: case C_SV:
282: return C_vs;
283: case C_SS:
284: return C_ss;
285: case C_VV:
286: return C_vv;
287: }
288: }
289:
290: // define width of a node
291: static inline int widthof(int node)
292: {
293: int width = max(Width[node],Nodesep);
294: if(N_self_edges > 0)
295: width += Selfwidth[node];
296: return width;
297: }
298:
299: // solve LP to get the positions
300: static void bestposition()
301: {
302: // calculate # of rows and cols
303: int jumpcnt = N_nodes - Maxlevel - 1;
304: int rows = N_edges + jumpcnt;
305: int cols = N_nodes + 2*N_edges + jumpcnt + 1;
306:
307: // allocate simplex arrays
308: long *weight = new long[N_edges];
309: long *objf = new long[cols];
310: long *auxf = new long[cols];
311: long **arrp = new long*[rows];
312: for(int i = 0; i < rows; ++i)
313: arrp[i] = new long[cols];
314:
315: // fill in first N_edges rows of the constraint array
316: int k = 0;
317: for(i = 0; i < N_nodes; i++)
318: {
319: for(edge_t *e = Outedges[i]; e; e = e->next)
320: {
321: arrp[k][N_nodes+k] = 1;
322: arrp[k][N_nodes+N_edges+k] = -1;
323: arrp[k][e->node] = -1;
324: arrp[k][i] = 1;
325: weight[k] = e->weight*edgefactor(i,e->node);
326: k++;
327: }
328: }
329:
330: // fill in last jumpcnt rows of constraint array
331: k = N_edges;
332: for(i = 0; i <= Maxlevel; i++)
333: {
334: for(int j = 0; j < N_level[i]-1; j++)
335: {
336: int v = Rank[i][j];
337: int w = Rank[i][j+1];
338: int vwidth = widthof(v);
339: int wwidth = widthof(w);
340:
341: arrp[k][cols-1] = (vwidth+wwidth)/2 + Nodesep;
342: arrp[k][v] = -1;
343: arrp[k][w] = 1;
344: arrp[k][N_nodes+N_edges+k] = -1;
345:
346: k++;
347: }
348: }
349:
350: // fill in last columns of objective function.
351: k = N_nodes+N_edges;
352: for(i = 0; i < N_edges; i++, k++)
353: objf[k] = 2*weight[i];
354:
355: // fill in last columns of auxiliary function
356: for(i = 0; i < jumpcnt; i++)
357: auxf[N_nodes + 2*N_edges + i] = 1;
358: auxf[cols-1] = -jumpcnt;
359:
360: // fill in first N_nodes columns of objective and aux functions
361: for(i = 0; i < N_nodes; i++)
362: {
363: long sum = 0;
364: for(k = 0; k < N_edges; ++k)
365: sum += arrp[k][i]*weight[k];
366: objf[i] = -sum;
367:
368: sum = 0;
369: for(int j = 0; j < jumpcnt; j++)
370: sum += arrp[N_edges + j][i];
371: auxf[i] = -sum;
372: }
373:
374: dag_simplex(rows,cols,arrp,auxf,objf,N_nodes,Nodepos);
375:
376: // restore space used by simplex
377: delete weight;
378: delete objf;
379: delete auxf;
380: for(i = 0; i < rows; ++i)
381: delete arrp[i];
382: delete arrp;
383: }
384:
385:
386: // try to straighten parts of a long edge
387: // the main for loop works from outside in to preserve any inherent symmetry
388: // in the part of the edge being worked on.
389:
390: static int straight(int n_nodes,int *nodes,int *minx,int *maxx,int *pos)
391: {
392: int n_done = 0;
393: if(n_nodes > 4)
394: {
395: int top = nodes[0], bot = nodes[n_nodes-1];
396: int top_lx = Nodepos[top] - Width[top]/2;
397: int top_rx = Nodepos[top] + Width[top]/2;
398: int bot_lx = Nodepos[bot] - Width[bot]/2;
399: int bot_rx = Nodepos[bot] + Width[bot]/2;
400: int d_top = Nodepos[top] - Nodepos[nodes[1]];
401: int d_bot = Nodepos[bot] - Nodepos[nodes[n_nodes-2]];
402: if(d_top*d_bot > 0)
403: {
404: for(int t = 1; t < n_nodes-1; ++t)
405: if(d_top*(Nodepos[nodes[t]]-Nodepos[nodes[t+1]]) <= 0)
406: break;
407: for(int b = n_nodes-2; b > t; --b)
408: if(d_bot*(Nodepos[nodes[b]]-Nodepos[nodes[b-1]]) <= 0)
409: break;
410: if(Nodepos[nodes[t]] >= top_lx && Nodepos[nodes[t]] <= top_rx)
411: t = 0;
412: if(Nodepos[nodes[b]] >= bot_lx && Nodepos[nodes[b]] <= bot_rx)
413: b = n_nodes-1;
414: if(t > 0 || b < n_nodes-1)
415: {
416: if(t >= 2)
417: n_done = straight(t+1,nodes,minx,maxx,pos);
418: if(b <= n_nodes-3)
419: n_done += straight(n_nodes-b,nodes+b,
420: minx+b,maxx+b,pos);
421: if((b-t) >= 2)
422: n_done += straight(b-t+1,nodes+t,
423: minx+t,maxx+t,pos);
424: return n_done;
425: }
426: }
427: }
428: // try to straighten all sub-intervals
429: for(int t = 0; t < n_nodes-2;)
430: {
431: for(int b = n_nodes-1; b > t+1; --b)
432: {
433: int vt = nodes[t];
434: int vb = nodes[b];
435: int xt = Nodepos[vt];
436: int xb = Nodepos[vb];
437: int yt = Levelpos[Level[vt]];
438: int yb = Levelpos[Level[vb]];
439:
440: // interpolate top/bottom nodes and compute new positions
441: int linear = 1;
442: for(int i = t+1, lev = Level[nodes[i]]; i < b; ++i, ++lev)
443: {
444: // new position of point
445: int cury = Levelpos[lev];
446: int curx = Nodepos[nodes[i]];
447: int newx = (xb-xt)*(cury-yb);
448: newx = (int)(((double)newx)/(yb-yt) + xb + .5);
449:
450: // easy case
451: if(newx == curx)
452: {
453: pos[i] = curx;
454: continue;
455: }
456:
457: // node width
458: int w = widthof(nodes[i])/2;
459:
460: // check constraints
461: int violated = 0;
462: if(newx-w < minx[i] || newx+w > maxx[i])
463: violated = 1;
464:
465: // make sure that yanking this straight
466: // won't cause node intersections
467: if(!violated)
468: {
469: int y = cury - Levelheight[lev]/2;
470: int x = (xb-xt)*(y-yb);
471: x = (int)(((double)x)/(yb-yt) + xb + .5);
472: if(x-w < minx[i] || x+w > maxx[i])
473: violated = 1;
474: }
475: if(!violated)
476: {
477: int y = cury + Levelheight[lev];
478: int x = (xb-xt)*(y-yb);
479: x = (int)(((double)x)/(yb-yt) + xb + .5);
480: if(x-w < minx[i] || x+w > maxx[i])
481: violated = 1;
482: }
483: if(violated && (b-t) > 2)
484: break;
485: else if(violated)
486: {
487: linear = 0;
488: int midx = (minx[i]+maxx[i])/2;
489: if((curx < midx && midx < newx) ||
490: (curx > midx && midx > newx))
491: newx = midx;
492: if(newx-w < minx[i] || newx+w > maxx[i])
493: break;
494: }
495: pos[i] = newx;
496: }
497:
498: // if the whole segment can be straightened
499: if(i == b)
500: {
501: for(i = t+1; i < b; ++i)
502: {
503: Nodepos[nodes[i]] = pos[i];
504: Deleted[nodes[i]] = linear ? 1 : 0;
505: }
506: Deleted[nodes[t]] = Deleted[nodes[b]] = 0;
507: n_done += b-t;
508:
509: // back to main loop
510: break;
511: }
512: }
513:
514: // start with the rest of the interval
515: t = b > t+1 ? b-1 : b;
516: }
517: return n_done;
518: }
519:
520:
521: // run thru all the edges once and try to straighten them out if
522: // no constraints are violated
523: static void setbounds(int *nodes, int n_nodes, int *minx, int *maxx)
524: {
525: for(int k = 1; k < n_nodes-1; ++k)
526: {
527: int v = nodes[k];
528: int i = Invert[v];
529: int *rank = Rank[Level[v]];
530: int a = rank[i+1];
531: maxx[k] = a < 0 ? Maxint : (Nodepos[a]-widthof(a)/2) - 5*Nodesep/4;
532: a = i > 0 ? rank[i-1] : -1;
533: minx[k] = a < 0 ? -Maxint : (Nodepos[a]+widthof(a)/2) + 5*Nodesep/4;
534: }
535: }
536:
537: static void smoothedges()
538: {
539: int *nodes = new int[Maxlevel+1];
540: int *pos = new int[Maxlevel+1];
541: int *minx = new int[Maxlevel+1];
542: int *maxx = new int[Maxlevel+1];
543:
544: int counter = 0;
545: int improving = 1;
546: while(improving)
547: {
548: improving = 0;
549: int dir = counter%2 ? -1 : 1;
550: for(int i = 0; i < Maxlevel; ++i)
551: {
552: int *begp = Rank[i];
553: int *endp = Rank[i]+N_level[i]-1;
554: int *vp = counter%2 ? endp : begp;
555: for(; vp >= begp && vp <= endp; vp += dir)
556: {
557: if(*vp >= N_real_nodes)
558: continue;
559: for(edge_t *e = Outedges[*vp]; e; e = e->next)
560: {
561: if(e->node < N_real_nodes)
562: continue;
563:
564: // get all the nodes
565: int n_nodes = 0;
566: nodes[n_nodes++] = *vp;
567: for(edge_t *f = e; f; f = f->chain)
568: nodes[n_nodes++] = f->node;
569:
570: // bounds for intermediate nodes
571: setbounds(nodes,n_nodes,minx,maxx);
572: int n = straight(n_nodes,nodes,minx,maxx,pos);
573: if(n > 2*n_nodes/3)
574: {
575: improving++;
576: e->node = -(e->node);
577: }
578: // need these to separate parallel edges
579: if(n > 0 && e->count > 1)
580: Deleted[nodes[1]] =
581: Deleted[nodes[n_nodes-2]] = 0;
582: }
583: }
584: }
585: }
586:
587: // reset signs of nodes
588: for(int v = 0; v < N_real_nodes; ++v)
589: for(edge_t *e = Outedges[v]; e; e = e->next)
590: if(e->node < 0)
591: e->node = -(e->node);
592: delete nodes;
593: delete pos;
594: delete minx;
595: delete maxx;
596: }
597:
598:
599: /*
600: Fast routines to assign positions.
601: The basic idea is to iteratively assign node positions based
602: on the median(s) of adjacent nodes. This heuristic is aided
603: with other heuristics that pack nodes closer based on local
604: improvements.
605: */
606:
607: // compute total edge length for a set of nodes
608: static int rlength(int *nodepos, int *node)
609: {
610: int length = 0;
611: for(; *node != -1; ++node)
612: for(edge_t *e = Outedges[*node]; e; e = e->next)
613: length += iabs(nodepos[*node]-nodepos[e->node])*e->weight;
614: if(length < 0)
615: panic("Integer overflow");
616: return length;
617: }
618:
619: // compute the length of all edges incident on a node givenits position
620: static int vlength(int v, int pos, int *nodepos)
621: {
622: int length = 0;
623: for(edge_t *e = Outedges[v]; e; e = e->next)
624: length += iabs(pos-nodepos[e->node])*e->weight;
625: for(e = Inedges[v]; e; e = e->next)
626: length += iabs(pos-nodepos[e->node])*e->weight;
627: return length;
628: }
629:
630:
631: // structures and functions for median computations
632: static struct MEDIAN
633: {
634: int position;
635: int weight;
636: int node;
637: };
638: static MEDIAN *Median;
639:
640: static struct PAIR
641: {
642: int left;
643: int right;
644: };
645:
646: static int *Priority;
647: static int prioritycmp(register int *one, register int *two)
648: {
649: return Priority[*two]-Priority[*one];
650: }
651:
652:
653: // compute the weighted median of a bunch of nodes
654: static int mediancmp(register MEDIAN *one, register MEDIAN *two)
655: {
656: return one->position - two->position;
657: }
658: static void wmedian(int n, PAIR *p)
659: {
660: if(n <= 0)
661: {
662: p->left = p->right = -1;
663: return;
664: }
665: if(n == 1)
666: {
667: p->left = p->right = 0;
668: return;
669: }
670:
671: // median weight;
672: int w_total = 0;
673: for(int i = 0; i < n; ++i)
674: w_total += Median[i].weight;
675: int w_median = (w_total+1)/2;
676:
677: // sort and the median
678: qsort((char*)Median,n,sizeof(Median[0]),(qsortcmpfun)mediancmp);
679: int w_accum = 0;
680: for(i = 0; i < n; ++i)
681: {
682: w_accum += Median[i].weight;
683: if(w_accum >= w_median)
684: break;
685: }
686: if(i == n)
687: panic("Bad median");
688: if(i < n-1 && (w_total%2) == 0 && w_accum == w_median)
689: {
690: p->left = i;
691: p->right = i+1;
692: }
693: else p->left = p->right = i;
694: }
695:
696: // assign position for nodes in a rank relative to an adjacent rank
697: static void assignpos(int *pos, int *node, int *rank, PAIR *median, int *priority)
698: {
699: for(; *node != -1; ++node)
700: {
701: // boundary in which the node must appear
702: int maxx = Maxint, minx = -Maxint;
703: int here = Invert[*node];
704: for(register int k = here+1; rank[k] != -1; ++k)
705: if(priority[rank[k]] > priority[*node])
706: break;
707: if(rank[k] != -1)
708: {
709: maxx = pos[rank[k]];
710: for(int n = k-1; n > here; --n)
711: maxx -= widthof(rank[n]);
712: maxx -= (k-here)*Nodesep;
713: maxx -= (widthof(*node)+widthof(rank[k]))/2;
714: }
715: for(k = here-1; k >= 0; --k)
716: if(priority[rank[k]] > priority[*node])
717: break;
718: if(k >= 0)
719: {
720: minx = pos[rank[k]];
721: for(int n = k+1; n < here; ++n)
722: minx += widthof(rank[n]);
723: minx += (here-k)*Nodesep;
724: minx += (widthof(*node)+widthof(rank[k]))/2;
725: }
726: if(maxx <= minx)
727: {
728: pos[*node] = minx;
729: continue;
730: }
731:
732: // compute new position
733: int set = 0;
734: PAIR *p = median + *node;
735: if(p->left >= 0)
736: {
737: int newx = (pos[p->left]+pos[p->right]+1)/2;
738: if(newx >= minx && newx <= maxx)
739: {
740: pos[*node] = newx;
741: ++set;
742: }
743: }
744: if(!set)
745: {
746: int curx = pos[*node];
747: if(curx < minx)
748: pos[*node] = minx;
749: else if(curx > maxx)
750: pos[*node] = maxx;
751: //else pos[*node] = curx;
752: }
753: }
754: }
755:
756:
757: // reassign positions using medians.
758: // Nodes are considered in decreasing priority
759: static void medianpos(int *pos, int **order, int *priority, PAIR *median, int isdown)
760: {
761: for(int l = 1; l <= Maxlevel; ++l)
762: {
763: int lev = isdown ? l : Maxlevel-l;
764: assignpos(pos,order[l],Rank[l],median,priority);
765: }
766: }
767:
768:
769: // same as medianpos, but for edges
770: static void minedge(int *nodepos)
771: {
772: for(int v = 0; v < N_real_nodes; ++v)
773: for(edge_t *e = Outedges[v]; e; e = e->next)
774: {
775: int w = e->node;
776: if(w >= N_real_nodes || nodepos[w] != nodepos[v])
777: continue;
778:
779: // interval in which the edge can move
780: int v_min, v_max;
781: int i = Invert[v];
782: int *rank = Rank[Level[v]];
783: int a = i > 0 ? rank[i-1] : -1;
784: v_min = a<0 ? -Maxint : nodepos[a]+(widthof(v)+widthof(a))/2 + Nodesep;
785: a = rank[i+1];
786: v_max = a<0 ? Maxint : nodepos[a]-(widthof(v)+widthof(a))/2 - Nodesep;
787: rank = Rank[Level[v]+1];
788:
789: int minx, maxx;
790: i = Invert[w];
791: a = i > 0 ? rank[i-1] : -1;
792: minx = a<0 ? -Maxint : nodepos[a]+(widthof(w)+widthof(a))/2 + Nodesep;
793: a = rank[i+1];
794: maxx = a<0 ? Maxint : nodepos[a]-(widthof(w)+widthof(a))/2 - Nodesep;
795: minx = max(minx,v_min);
796: maxx = min(maxx,v_max);
797: if(maxx <= minx)
798: continue;
799:
800: // make list of adjacent nodes
801: MEDIAN *m = Median;
802: int count = 0;
803: for(int top = 0; top < 2; ++top)
804: {
805: int here = top ? v : w;
806: int there = top ? w : v;
807: for(int in = 0; in < 2; ++in)
808: {
809: edge_t *e = in ? Inedges[here] : Outedges[here];
810: for(; e; e = e->next)
811: {
812: if(e->node == there)
813: continue;
814: m->position = nodepos[e->node];
815: m->weight = e->weight;
816: ++m;
817: ++count;
818: }
819: }
820: }
821: if(count <= 0)
822: continue;
823:
824: // get the median(s)
825: PAIR p;
826: wmedian(count,&p);
827: int newx = (Median[p.left].position+Median[p.right].position+1)/2;
828: if(newx >= minx && newx <= maxx)
829: {
830: nodepos[v] = newx;
831: nodepos[w] = newx;
832: }
833: }
834: }
835:
836:
837: // local improvement for each node based on its entire adjacent nodes
838: // This is the same as doing a descent from the current assignment to
839: // improve the objective function. This is why the algorithm will
840: // terminate.
841: static void minnode(int *nodepos, int *marked, int *queue)
842: {
843: // local inline functions
844: #define ENQUEUE(x) queue[tail] = x; if(++tail >= N_nodes+1) tail = 0;
845: #define DEQUEUE(x) x = queue[head]; if(++head >= N_nodes+1) head = 0;
846: #define NOTNULL() head != tail
847:
848: // marked[] and queue[] have the set of nodes to be considered
849: for(int n = 0; n < N_nodes; ++n)
850: {
851: marked[n] = 1;
852: queue[n] = n;
853: }
854:
855: int head = 0, tail = N_nodes;
856: while(NOTNULL())
857: {
858: int here;
859:
860: DEQUEUE(here);
861: marked[here] = 0;
862:
863: // boundary in which the node must appear
864: int here_w = widthof(here);
865: n = Invert[here];
866: int *rank = Rank[Level[here]];
867: int l = n <= 0 ? -1 : rank[n-1];
868: int r = rank[n+1];
869: int minx = l < 0 ? -Maxint : nodepos[l]+Nodesep+(here_w+widthof(l))/2;
870: int maxx = r < 0 ? Maxint : nodepos[r]-Nodesep-(here_w+widthof(r))/2;
871: if(maxx <= minx)
872: continue;
873:
874: int cnt = 0;
875: MEDIAN *median = Median;
876: for(int down = 0; down < 2; ++down)
877: {
878: edge_t *e = down ? Outedges[here] : Inedges[here];
879: for(; e; e = e->next)
880: {
881: median->position = nodepos[e->node];
882: median->weight = e->weight;
883: median->node = -1;
884: cnt++;
885: median++;
886: }
887: }
888: if(cnt == 0)
889: continue;
890:
891: PAIR p;
892: wmedian(cnt,&p);
893: int newx = (Median[p.left].position+Median[p.right].position+1)/2;
894: int correction = 0;
895: if(newx < minx || newx > maxx)
896: {
897: correction = 1;
898: newx = newx < minx ? minx : maxx;
899: }
900: if(newx == nodepos[here])
901: continue;
902:
903: int newlength = vlength(here,newx,nodepos);
904: int oldlength = vlength(here,nodepos[here],nodepos);
905: if(newlength > oldlength || (newlength == oldlength && correction))
906: continue;
907:
908: // update position and mark neighbors for consideration
909: nodepos[here] = newx;
910: if(newlength < oldlength)
911: {
912: for(down = 0; down < 2; ++down)
913: {
914: int adj = down ? l : r;
915: if(adj >= 0 && !marked[adj])
916: {
917: marked[adj] = 1;
918: ENQUEUE(adj);
919: }
920: edge_t *e = down ? Outedges[here] : Inedges[here];
921: for(; e; e = e->next)
922: {
923: adj = e->node;
924: if(marked[adj])
925: continue;
926: marked[adj] = 1;
927: ENQUEUE(adj);
928: }
929: }
930: }
931: }
932: #undef ENQUEUE
933: #undef DEQUEUE
934: #undef NOTNULL
935: }
936:
937: // pack cuts that are not tight
938: static int *_nodepos;
939: static int rightcmp(int *p1, int *p2)
940: {
941: int right1 = _nodepos[*p1] + widthof(*p1)/2;
942: int right2 = _nodepos[*p2] + widthof(*p2)/2;
943: return right1-right2;
944: }
945: static int leftcmp(int *p1, int *p2)
946: {
947: int left1 = _nodepos[*p1] - widthof(*p1)/2;
948: int left2 = _nodepos[*p2] - widthof(*p2)/2;
949: return left1-left2;
950: }
951: static void verticalcut(int *nodepos, int *nodes, int maxwidth)
952: {
953: // sort vertices by their positions
954: _nodepos = nodepos;
955: qsort((char*)nodes,N_nodes,sizeof(nodes[0]),(qsortcmpfun)rightcmp);
956:
957: // go thru each possible cut
958: for(int i = 0; i+1 < N_nodes; ++i)
959: {
960: int v = nodes[i];
961: int cut = nodepos[v] + widthof(v)/2 + Nodesep;
962:
963: int leftside = Maxint;
964: for(int k = i+1; k < N_nodes; ++k)
965: {
966: int w = nodes[k];
967:
968: // don't need to check further
969: if((nodepos[w]-2*(maxwidth+Nodesep)) > leftside)
970: break;
971:
972: // define leftside
973: int wleft = nodepos[w] - widthof(w)/2;
974: if(wleft < leftside)
975: leftside = wleft;
976: if(leftside <= cut)
977: break;
978: }
979: // can improve
980: int shift = leftside-cut;
981: if(shift > 0)
982: for(k = i+1; k < N_nodes; ++k)
983: nodepos[nodes[k]] -= shift;
984: }
985: }
986:
987: static void mincut(int *nodepos, int *marked, int *queue, int *nodes, int dir)
988: {
989: #define ENQUEUE(x) queue[tail] = x; if(++tail >= N_nodes) tail = 0;
990: #define DEQUEUE(x) x = queue[head]; if(++head >= N_nodes) head = 0;
991: #define NOTNULL() head != tail
992:
993: _nodepos = nodepos;
994: qsort((char*)nodes,N_nodes,sizeof(nodes[0]),(dir < 0 ? (qsortcmpfun)leftcmp : (qsortcmpfun)rightcmp));
995: int cut, lastcut = Maxint;
996: int n = dir > 0 ? 0 : N_nodes-1;
997: for(; n >= 0 && n < N_nodes; n += dir, lastcut = cut)
998: {
999: if(nodes[n] >= N_real_nodes)
1000: continue;
1001:
1002: // where the graph is being cut
1003: cut = nodepos[nodes[n]]+dir*(widthof(nodes[n])/2+Nodesep);
1004: if((dir > 0 && cut <= lastcut) || (dir < 0 && cut >= lastcut))
1005: continue;
1006:
1007: // try to shift connected components in the right of the cut
1008: for(int i = 0; i < N_nodes; ++i)
1009: marked[i] = 0;
1010: int component = 1;
1011: for(int m = n+dir; m >= 0 && m < N_nodes; m += dir)
1012: {
1013: int v = nodes[m];
1014: if(marked[v])
1015: continue;
1016: int side = nodepos[v] - dir*widthof(v)/2;
1017: if((dir > 0 && side <= cut) || (dir < 0 && side >= cut))
1018: continue;
1019:
1020: // bf-search for all nodes in this component
1021: int head = 0, tail = 0;
1022: ENQUEUE(v);
1023: marked[v] = component;
1024: while(NOTNULL())
1025: {
1026: int here;
1027: DEQUEUE(here);
1028: for(int in = 0; in < 2; ++in)
1029: {
1030: edge_t *e = in ? Inedges[here]:Outedges[here];
1031: for(; e; e = e->next)
1032: {
1033: int w = e->node;
1034: if(marked[w])
1035: continue;
1036: int side = nodepos[w]-dir*widthof(w)/2;
1037: if((dir > 0 && side <= cut) ||
1038: (dir < 0 && side >= cut))
1039: continue;
1040: ENQUEUE(w);
1041: marked[w] = component;
1042: }
1043: }
1044: }
1045:
1046: int shift = Maxint;
1047: for(i = 0; i < tail; ++i)
1048: {
1049: int v = queue[i];
1050: int a = Invert[v];
1051: if(dir > 0)
1052: a = a <= 0 ? -1 : Rank[Level[v]][a-1];
1053: else a = Rank[Level[v]][a+1];
1054: if(a < 0 || marked[v] == marked[a])
1055: continue;
1056: int side = nodepos[a] + dir*(widthof(a)/2+Nodesep);
1057: if((dir > 0 && side < cut) || (dir < 0 && side > cut))
1058: side = cut;
1059: int slack = nodepos[v] - dir*(widthof(v)/2);
1060: slack = dir > 0 ? slack-side : side-slack;
1061: if(slack < shift)
1062: shift = slack;
1063: if(shift <= 0)
1064: break;
1065: }
1066: if(shift < Maxint && shift > 0)
1067: for(i = 0; i < tail; ++i)
1068: nodepos[queue[i]] -= dir*shift;
1069: component++;
1070: }
1071: }
1072: }
1073:
1074: // straighten paths
1075: static void getminmax(int *nodepos, int v, int& minx, int& maxx)
1076: {
1077: int *rank = Rank[Level[v]];
1078: int i = Invert[v];
1079: int l = i > 0 ? rank[i-1] : -1;
1080: int r = rank[i+1];
1081: minx = l < 0 ? -Maxint : nodepos[l] + (widthof(l)+widthof(v))/2 + Nodesep;
1082: maxx = r < 0 ? Maxint : nodepos[r] - (widthof(r)+widthof(v))/2 - Nodesep;
1083: }
1084:
1085: static void minpath(int *nodepos, int *nodes)
1086: {
1087: for(int lev = 0; lev < Maxlevel; ++lev)
1088: for(int *rank = Rank[lev]; *rank >= 0; ++rank)
1089: for(edge_t *e = Outedges[*rank]; e; e = e->next)
1090: {
1091: // find a straightenable path from v
1092: int n_nodes = 0;
1093: int minx = -Maxint, maxx = Maxint;
1094: int v = *rank, w = e->node;
1095: if(pathdegree(v) == 1 && !Inedges[v])
1096: {
1097: nodes[n_nodes++] = v;
1098: getminmax(nodepos,v,minx,maxx);
1099: }
1100: edge_t *f = e;
1101: while(1)
1102: {
1103: w = f->node;
1104: if(pathdegree(w) > 2)
1105: break;
1106:
1107: int lx, rx;
1108: getminmax(nodepos,w,lx,rx);
1109: if(lx > maxx || rx < minx)
1110: break;
1111:
1112: minx = max(minx,lx);
1113: maxx = min(maxx,rx);
1114: nodes[n_nodes++] = w;
1115: if(!Outedges[w])
1116: break;
1117: else f = Outedges[w];
1118: }
1119: if(n_nodes <= 0 || minx > maxx)
1120: continue;
1121:
1122: // new position
1123: int vx = nodepos[v];
1124: int wx = nodepos[w];
1125: vx = vx < minx ? minx : (vx > maxx ? maxx : vx);
1126: wx = wx < minx ? minx : (wx > maxx ? maxx : wx);
1127: int newx;
1128: if(v == nodes[0] && w != nodes[n_nodes-1])
1129: newx = wx;
1130: else if(w == nodes[n_nodes-1] && v != nodes[0])
1131: newx = vx;
1132: else newx = e->weight > f->weight ? vx :
1133: (f->weight > e->weight ? wx : (vx+wx+1)/2);
1134: for(int i = 0; i < n_nodes; ++i)
1135: nodepos[nodes[i]] = newx;
1136: }
1137: }
1138:
1139:
1140: // define median pointers to help medianpos().
1141: // The observation here is that since the relative order of
1142: // adjacent nodes of a node in the next or last level is always
1143: // maintained, the median node is known regardless of its position.
1144: // So we compute it once and avoid doing qsorts during position
1145: // reassignment.
1146:
1147: static void define_median(PAIR *median, edge_t **edges)
1148: {
1149: for(int i = 0; i < N_nodes; ++i)
1150: {
1151: int cnt = 0;
1152: MEDIAN *mp = Median;
1153: for(edge_t *e = edges[i]; e; e = e->next)
1154: {
1155: mp->weight = e->weight;
1156: mp->position = Invert[e->node];
1157: mp->node = e->node;
1158: ++mp;
1159: ++cnt;
1160: }
1161: if(cnt == 0)
1162: median[i].left = median[i].right = -1;
1163: else
1164: {
1165: PAIR p;
1166: wmedian(cnt,&p);
1167: median[i].left = Median[p.left].node;
1168: median[i].right = Median[p.right].node;
1169: }
1170: }
1171: }
1172:
1173: // compute starting positions
1174: static void startpos(int *pos, int **u_order, PAIR *u_median, int *u_prior,
1175: int **d_order, PAIR *d_median, int *d_prior)
1176: {
1177: int maxrank = -1, maxwidth = 0;
1178: for(int i = 0; i <= Maxlevel; ++i)
1179: {
1180: int wlast = 0;
1181: int xlast = -Nodesep;
1182: for(int *rank = Rank[i]; *rank >= 0; ++rank)
1183: {
1184: int wnow = widthof(*rank);
1185: xlast += Nodesep+(wlast+wnow)/2;
1186: pos[*rank] = xlast;
1187: wlast = wnow;
1188: }
1189: if(xlast+wlast/2 > maxwidth)
1190: {
1191: maxrank = i;
1192: maxwidth = xlast+wlast/2;
1193: }
1194: }
1195: for(i = maxrank-1; i >= 0; --i)
1196: assignpos(pos,u_order[i],Rank[i],u_median,u_prior);
1197: for(i = maxrank+1; i <= Maxlevel; ++i)
1198: assignpos(pos,d_order[i],Rank[i],d_median,d_prior);
1199: }
1200:
1201: // define orders
1202: static void define_order(int *u_prior, int **u_order, int *d_prior, int **d_order)
1203: {
1204: for(int i = 0; i <= Maxlevel; ++i)
1205: {
1206: int n_rank = N_level[i];
1207: int *rank = Rank[i];
1208: int *up = u_order[i];
1209: int *down = d_order[i];
1210:
1211: for(int k = n_rank; k >= 0; --k)
1212: up[k] = down[k] = rank[k];
1213: for(k = 0; rank[k] != -1;)
1214: {
1215: for(int ek = k+1; rank[ek] != -1; ++ek)
1216: if(Component[rank[ek]] != Component[rank[ek-1]])
1217: break;
1218: Priority = u_prior;
1219: qsort((char*)(up+k),ek-k,sizeof(up[0]),(qsortcmpfun)prioritycmp);
1220: Priority = d_prior;
1221: qsort((char*)(down+k),ek-k,sizeof(down[0]),(qsortcmpfun)prioritycmp);
1222:
1223: k = ek;
1224: }
1225: for(k = 0; rank[k] != -1; ++k, --n_rank)
1226: {
1227: u_prior[up[k]] = n_rank;
1228: d_prior[down[k]] = n_rank;
1229: }
1230: }
1231: }
1232:
1233: static void goodposition()
1234: {
1235: const int Maxiter = 8, Minquit = 4;
1236:
1237: // set priorities of nodes
1238: int *d_prior = new int[N_nodes];
1239: int *u_prior = new int[N_nodes];
1240: for(int i = 0; i < N_nodes; ++i)
1241: {
1242: d_prior[i] = u_prior[i] = 0;
1243: for(edge_t *e = Inedges[i]; e; e = e->next)
1244: d_prior[i] += (e->weight *= edgefactor(e->node,i));
1245: for(e = Outedges[i]; e; e = e->next)
1246: u_prior[i] += (e->weight *= edgefactor(i,e->node));
1247: }
1248:
1249: // construct the order in which positions are assigned to nodes
1250: int **u_order = new int*[Maxlevel+1];
1251: int **d_order = new int*[Maxlevel+1];
1252: for(i = 0; i <= Maxlevel; ++i)
1253: {
1254: u_order[i] = new int[N_level[i]+1];
1255: d_order[i] = new int[N_level[i]+1];
1256: }
1257: define_order(u_prior,u_order,d_prior,d_order);
1258:
1259: // compute the median pointers
1260: Median = new MEDIAN[N_nodes];
1261: int *nodepos = new int[N_nodes];
1262: PAIR *u_median = new PAIR[N_nodes];
1263: PAIR *d_median = new PAIR[N_nodes];
1264: define_median(d_median,Inedges);
1265: define_median(u_median,Outedges);
1266:
1267: // temp space for graph search
1268: int *marked = new int[N_nodes];
1269: int *queue = new int[N_nodes+1];
1270: int *nodes = new int[N_nodes];
1271: int maxwidth = 0;
1272: for(i = 0; i < N_nodes; ++i)
1273: {
1274: nodes[i] = i;
1275: maxwidth = max(maxwidth,widthof(i));
1276: }
1277:
1278: long bestlength = -1;
1279: int improved = 0;
1280: for(int upfirst = 1; upfirst >= 0; --upfirst)
1281: {
1282: // compute starting positions
1283: startpos(nodepos,u_order,u_median,u_prior,d_order,d_median,d_prior);
1284: int counter = 1;
1285: for(int iter = 0; iter < Maxiter; ++iter)
1286: {
1287: if(upfirst)
1288: {
1289: if(iter%2 == 0)
1290: medianpos(nodepos,u_order,u_prior,u_median,0);
1291: else medianpos(nodepos,d_order,d_prior,d_median,1);
1292: }
1293: else
1294: {
1295: if(iter%2 == 0)
1296: medianpos(nodepos,d_order,d_prior,d_median,1);
1297: else medianpos(nodepos,u_order,u_prior,u_median,0);
1298: }
1299:
1300: // local optimization
1301: minedge(nodepos);
1302: minnode(nodepos,marked,queue);
1303: minpath(nodepos,queue);
1304: verticalcut(nodepos,nodes,maxwidth);
1305:
1306: long curlength = 0;
1307: for(i = 0; i < Maxlevel; ++i)
1308: curlength += rlength(nodepos,Rank[i]);
1309: if(bestlength < 0 || curlength < bestlength)
1310: {
1311: bestlength = curlength;
1312: for(i = 0; i < N_nodes; ++i)
1313: Nodepos[i] = nodepos[i];
1314: improved++;
1315: counter = 1;
1316: }
1317: else if(counter++ > Minquit)
1318: break;
1319: }
1320: }
1321:
1322: mincut(Nodepos,marked,queue,nodes,1);
1323: mincut(Nodepos,marked,queue,nodes,-1);
1324: minedge(Nodepos);
1325: minnode(Nodepos,marked,queue);
1326: minpath(Nodepos,queue);
1327:
1328: // reclaim space
1329: for(i = 0; i <= Maxlevel; ++i)
1330: {
1331: delete u_order[i];
1332: delete d_order[i];
1333: }
1334: delete u_order;
1335: delete d_order;
1336: delete u_prior;
1337: delete d_prior;
1338: delete u_median;
1339: delete d_median;
1340: delete marked;
1341: delete queue;
1342: delete nodes;
1343: delete nodepos;
1344: delete Median;
1345: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.