|
|
researchv10 Norman
#include "draw_dag.h"
#include <sys/types.h>
#include <sys/times.h>
#define TIC 60
/*
Create an ordering of nodes on each level.
The ordering attempts to minimize edge crossing
*/
static float Convergence = .005; // convergence rate
static int Maxiter = 24; // max number of iterations
static int Minquit = 10; // quit after this many iterations
// if not improved
static int **List; // temp space for rank construction
static int *Count, // to aid rcross()
*Cross, // current numbers of crossings between levels
*Levelmarked; // used in transpose()
static edge_t **Flatedges; // transitive closure of flat edges
static int *Flatindeg; // in-degree of nodes using Flatedges
static long *Median;
static int *Medadj;
static void decompose(),
flatclosure(),
removestem(),
installstem(int**,int);
static int cross_stat(int**), rcross(int*,int*),
buildrank(int**,int**,int);
extern void inversion(int**);
void dag_ranks()
{
struct tms begtm, endtm;
int startcross = 0;
if(Verbose)
{
#ifndef TIMING
fprintf(stderr,"Minimize edge crossing\n");
#endif
times(&begtm);
}
#ifdef TIMING
Maxiter = 20;
Minquit = 20;
Convergence = .0001;
#else
// reset parameters for speed
if(N_real_nodes > 500)
{
float density = ((float) N_edges)/N_nodes;
if(density >= 2.)
{
Maxiter = N_real_nodes <= 1000 ? 8 : 4;
Minquit = N_real_nodes <= 1000 ? 4 : 2;
Convergence = N_real_nodes <= 1000 ? .01 : .1;
}
}
#endif
// number of nodes in each level
N_level = new int[Maxlevel+1];
for(int i = 0; i < N_nodes; i++)
N_level[Level[i]] += 1;
// create space for node lists
List = new int*[Maxlevel+1];
Rank = new int*[Maxlevel+1];
int **tryout = new int*[Maxlevel+1];
int **target = new int*[Maxlevel+1];
for(i = 0; i <= Maxlevel; ++i)
{
List[i] = new int[N_level[i]+1];
tryout[i] = new int[N_level[i]+1];
Rank[i] = new int[N_level[i]+1];
Rank[i][0] = -1;
}
// space for intermediate calculations
Invert = new int[N_nodes];
Count = new int[N_nodes];
Cross = new int[Maxlevel];
Levelmarked = new int[Maxlevel+1];
Median = new long[N_nodes];
Medadj = new int[N_nodes+N_repeat_edges];
// compute the transitive closure of the set of flat edges
if(N_flat_edges > 0)
{
Flatedges = new edge_t*[N_real_nodes];
Flatindeg = new int[N_real_nodes];
flatclosure();
}
// detect connected components
Component = new int[N_nodes];
decompose();
// build up the ordering component by component
for(int component = 1; component <= N_component; ++component)
{
for(i = 0; i <= Maxlevel; ++i)
{
for(int k = 0; ; ++k)
if(Rank[i][k] == -1)
break;
target[i] = Rank[i]+k;
}
startcross += buildrank(target,tryout,component);
}
inversion(Rank);
if(Verbose)
{
times(&endtm);
int total = (endtm.tms_utime - begtm.tms_utime)+
(endtm.tms_stime - begtm.tms_stime);
int percent = (total - (total/TIC)*TIC)*100/TIC;
int bestcross = cross_stat(Rank);
#ifdef TIMING
fprintf(stderr,"%d\t%d\t%d.%02d\t",
startcross,bestcross,total/TIC,percent);
#else
fprintf(stderr,"Starting crossing number = %d\n",startcross);
fprintf(stderr,"Crossings between ranks:\n");
for(i = 0; i < Maxlevel; ++i)
fprintf(stderr,"(%2d,%2d):\tx=%d\n",i,i+1,Cross[i]);
fprintf(stderr,"Ending crossing number = %d\n",bestcross);
fprintf(stderr,"Time in minimizing crossing: %d.%02ds\n",
total/TIC, percent);
#endif
}
// release temp space
delete Median;
delete Medadj;
delete Levelmarked;
delete Cross;
delete Count;
for(i = 0; i <= Maxlevel; ++i)
{
delete List[i];
delete tryout[i];
}
delete List;
delete tryout;
delete target;
if(N_flat_edges > 0)
deledges(N_real_nodes,Flatedges);
}
// compute the transitive closure of flat edges
// and remove resulting bidirectional edges.
// The final graph defines a partial order on the nodes.
static void flatedges(int root, int here, int *marked)
{
marked[here] = 1;
for(edge_t *e = Noedges[here]; e; e = e->next)
if(e->weight > 0 && !marked[e->node])
{
// check for bi-directional edges
edge_t *lastf = NULL, *f = Flatedges[e->node];
for(; f != NULL; lastf = f, f = f->next)
if(f->node == root)
break;
// remove bidirectional edges
if(f != NULL)
{
if(lastf)
lastf->next = f->next;
else Flatedges[e->node] = f->next;
delete f;
}
// add a new edge
else Flatedges[root] =
new edge_t(e->node,0,0,0,Flatedges[root]);
flatedges(root,e->node,marked);
}
}
static void flatclosure()
{
int *marked = Count;
for(int v = 0; v < N_real_nodes; ++v)
{
for(int i = 0; i < N_real_nodes; ++i)
marked[i] = 0;
flatedges(v,v,marked);
}
for(v = 0; v < N_real_nodes; ++v)
for(edge_t *e = Flatedges[v]; e; e = e->next)
Flatindeg[e->node]++;
}
// assign component numbers to nodes so that nodes in the same
// connected component get the same number.
static void decompose()
{
// make a reverse flat edge list
edge_t **noedges = NULL;
if(N_flat_edges > 0)
{
noedges = new edge_t*[N_real_nodes];
for(int v = 0; v < N_real_nodes; ++v)
for(edge_t *e = Noedges[v]; e; e = e->next)
{
int w = e->node;
if(w != v)
noedges[w] = new edge_t(v,0,0,0,noedges[w]);
}
}
// all adjacent edge lists
edge_t **adjlist[5];
adjlist[0] = Outedges;
adjlist[1] = Inedges;
adjlist[2] = N_flat_edges > 0 ? Noedges : NULL;
adjlist[3] = N_flat_edges > 0 ? noedges : NULL;
adjlist[4] = NULL;
// space for the breadth-first search
int *queue = new int[N_nodes];
N_component = 0;
for(int n = 0; n < N_nodes; ++n)
{
if(Component[n] > 0)
continue;
// initialize for component search
N_component += 1;
Component[n] = N_component;
// breadth-first search
int first = 0, last = 1;
queue[0] = n;
while(first < last)
{
int v = queue[first++];
for(int i = 0; adjlist[i] != NULL; ++i)
{
if(v >= N_real_nodes && i >= 2)
continue;
for(edge_t *e = adjlist[i][v]; e; e = e->next)
{
int w = e->node;
if(Component[w] != 0)
continue;
Component[w] = N_component;
queue[last++] = w;
}
}
}
}
delete queue;
if(N_flat_edges > 0)
deledges(N_real_nodes,noedges);
}
// copy a ranking to another
static void copyrank(int **from, int **to)
{
for(int i = 0; i <= Maxlevel; ++i)
{
register int *frp = from[i], *top = to[i];
while(*frp != -1)
*top++ = *frp++;
*top = -1;
}
}
// see if one vertex must be left of another
static int left2right(int v, int w)
{
if(v >= N_real_nodes || w >= N_real_nodes)
return 0;
for(edge_t *e = Flatedges[v]; e; e = e->next)
if(e->node == w)
return 1;
return 0;
}
// construct nodes reachable from 'here' in post-order.
// This is the same as doing a topological sort in reverse order.
static void postorder(int here,int *marked,int *list, int &n_search)
{
marked[here] = 1;
for(edge_t *e = Flatedges[here]; e; e = e->next)
if(!marked[e->node])
postorder(e->node,marked,list,n_search);
list[n_search++] = here;
}
static void reorder(int *list, int *marked, int *temp)
{
for(int i = 0; list[i] >= 0; ++i)
marked[list[i]] = 0;
int n_temp = 0;
for(i = 0; list[i] >= 0; ++i)
{
// dummy nodes
if(list[i] >= N_real_nodes)
temp[n_temp++] = list[i];
// real nodes
else if(!marked[list[i]] && Flatindeg[list[i]] == 0)
{
int n_search = 0;
postorder(list[i],marked,temp+n_temp,n_search);
int *left = temp+n_temp, *right = left+n_search-1;
for(; left < right; ++left, --right)
{
int t = *left;
*left = *right;
*right = t;
}
n_temp += n_search;
}
}
// recopy the list
for(i = 0; list[i] >= 0; ++i)
list[i] = temp[i];
}
// Compute an initial ordering of nodes belonging to some connected component.
static void initrank(int *rank[], int component, int isdown)
{
// use breadth-first search to build the list
int *marked = new int[N_nodes];
int *queue = new int[N_nodes];
int *n_lev = new int[Maxlevel+1];
// now do the rank partitioning by components
edge_t **those = isdown ? Inedges : Outedges;
edge_t **these = isdown ? Outedges : Inedges;
for(int n = 0; n < N_nodes; ++n)
{
if(marked[n] || those[n] || Component[n] != component ||
(n < N_real_nodes && N_flat_edges <= 0 && Stems[n]))
continue;
marked[n] = 1;
queue[0] = n;
int first = 0, last = 1;
while(first < last)
{
int node = queue[first++];
rank[Level[node]][n_lev[Level[node]]++] = node;
for(edge_t *e = these[node]; e; e = e->next)
{
if(!marked[e->node])
{
marked[e->node] = 1;
queue[last++] = e->node;
}
}
}
}
for(n = 0; n <= Maxlevel; ++n)
{
rank[n][n_lev[n]] = -1;
// make sure that implicit left-to-right orders are satisfied
if(N_flat_edges > 0)
reorder(rank[n],marked,queue);
}
delete marked;
delete queue;
delete n_lev;
}
// Compute the crossing of two levels.
// We'll assume that edges point from top to bot.
static int rcross(int *top, int *bot)
{
// clear up work space
for (int i = 0; bot[i] != -1; i++)
Count[i] = 0;
register int cross = 0; // the crossing count
int max = 0; // index we have to sum to during counting
for (i = 0; top[i] != -1; i++)
{
register edge_t *e;
if (max > 0)
{
/* iterate over each edge */
for (e = Outedges[top[i]]; e; e = e->next)
for(int k = Invert[e->node]+1; k <= max; k++)
cross += Count[k]*e->count;
}
/* update the counts */
for (e = Outedges[top[i]]; e; e = e->next)
{
register int inv = Invert[e->node];
if (inv > max)
max = inv;
Count[inv] += e->count;
}
}
return cross;
}
// compute the inverted indices of all vertices in their rank lists
// and the crossing statistics between ranks.
static void inversion(int **ranks)
{
for(int i = 0; i <= Maxlevel; ++i)
{
int *rank = ranks[i];
for(int k = 0; rank[k] != -1; ++k)
Invert[rank[k]] = k;
}
}
static int cross_stat(int **ranks)
{
int cross = 0;
for(int i = 0; i < Maxlevel; ++i)
cross += (Cross[i] = rcross(ranks[i],ranks[i+1]));
return cross;
}
// Compute the crossing of incident edges of two vertices
static int vcross(edge_t *adj1, edge_t *adj2)
{
register int cross = 0;
for (register edge_t *e2 = adj2; e2; e2 = e2->next)
{
register int inv = Invert[e2->node];
register int cnt = e2->count;
for(register edge_t *e1 = adj1; e1; e1 = e1->next)
if(Invert[e1->node] > inv)
cross += e1->count * cnt;
}
return cross;
}
// Reduce crossings by adjacent tranpositions of vertices
static void transpose(int **ranks, int revert)
{
// current number of crossings
int curcross = 0;
for(int l = 0; l < Maxlevel; ++l)
{
curcross += Cross[l];
Levelmarked[l] = 1;
}
if(curcross == 0)
return;
else Levelmarked[l] = 1;
while(1)
{
int delta = 0;
for (int lev = 0; lev <= Maxlevel; lev++)
if(Levelmarked[lev])
{
int *rank = ranks[lev];
if(rank[0] == -1 || rank[1] == -1)
{
Levelmarked[lev] = 0;
continue;
}
for(; rank[1] != -1; ++rank)
{
register int c0 = 0, c1 = 0;
int v = rank[0], w = rank[1];
if(N_flat_edges > 0 && left2right(v,w) != 0)
continue;
if(lev > 0)
{
edge_t *vadj = Inedges[v], *wadj = Inedges[w];
c0 += vcross(vadj,wadj);
c1 += vcross(wadj,vadj);
}
if(lev < Maxlevel)
{
edge_t *vadj = Outedges[v], *wadj = Outedges[w];
c0 += vcross(vadj,wadj);
c1 += vcross(wadj,vadj);
}
if (c1 < c0 || (c1 > 0 && c1 == c0 && revert))
{
Levelmarked[lev] = 2;
Invert[v]++;
Invert[w]--;
rank[0] = w;
rank[1] = v;
delta += c0-c1;
}
}
if(Levelmarked[lev] == 2)
{
Levelmarked[lev] = 1;
if(lev > 0)
{
Levelmarked[lev-1] = 1;
Cross[lev-1] = -1;
}
if(lev < Maxlevel)
{
Levelmarked[lev+1] = 1;
Cross[lev] = -1;
}
}
else Levelmarked[lev] = 0;
}
curcross -= delta;
float reduce = curcross <= 0 ? 0. : ((float)delta)/curcross;
if(curcross == 0 || reduce <= Convergence)
break;
}
// restore crossing statistics
for(l = 0; l < Maxlevel; ++l)
if(Cross[l] < 0)
Cross[l] = rcross(ranks[l],ranks[l+1]);
}
// sort nodes by their median values.
// This is basically a bubble sort that respects fixed points
static int mediansort(int *array, int nelt, int revert, int hasfixed)
{
int changed = 0;
register int *ep = array+nelt, *lp, *rp;
for(nelt -= 1; nelt > 0; --nelt)
{
lp = array;
while(lp < ep)
{
// -1 nodes are FIXED points
if(Median[*lp] < 0)
for(; lp < ep; lp++)
if(Median[*lp] >= 0)
break;
if(lp >= ep)
break;
// find the one that can be compared with
int muststay = 0;
for(rp = lp+1; rp < ep; ++rp)
{
if(N_flat_edges > 0 && left2right(*lp,*rp))
{
muststay = 1;
break;
}
if(Median[*rp] >= 0)
break;
}
if(rp >= ep)
break;
// switch if appropriate
if(!muststay &&
(Median[*lp] > Median[*rp] ||
(Median[*lp] == Median[*rp] && revert)))
{
register int t = *lp;
*lp = *rp;
*rp = t;
changed++;
}
lp = rp;
}
if(!hasfixed && !revert)
--ep;
}
return changed;
}
// reduce crossing using weighted median sorting
static int intcmp(register int *i1, register int *i2)
{
return *i1 - *i2;
}
static void median(int **ranks, int revert, int isup)
{
#define SCALE 16
for (int lev = 1; lev <= Maxlevel; lev++)
{
int thislev = isup ? Maxlevel-lev : lev;
int thatlev = isup ? thislev+1 : thislev-1;
edge_t **edges = isup ? Outedges : Inedges;
int *thislist = ranks[thislev];
int *thatlist = ranks[thatlev];
// compute medians of thislist
int hasfixed = 0;
for(int n_this = 0; thislist[n_this] != -1; n_this++)
{
register int node = thislist[n_this];
register int *endadj = Medadj;
for(register edge_t *e = edges[node]; e; e = e->next)
for(int k = e->count; k > 0; --k)
*endadj++ = Invert[e->node];
register int cnt = endadj-Medadj;
#ifdef BCENTER
// barycenter method
if(cnt > 0)
{
int ave = 0;
for(int m = 0; m < cnt; ++m)
ave += Medadj[m];
Median[node] = (SCALE*ave)/cnt;
}
else
{
Median[node] = -1;
hasfixed = 1;
}
#else /* a right median method */
#ifdef RMEDIAN
if(cnt > 1)
qsort((char*)Medadj,cnt,sizeof(Medadj[0]),(qsortcmpfun)intcmp);
if(cnt == 0)
{
Median[node] = -1;
hasfixed = 1;
}
else Median[node] = SCALE*Medadj[cnt/2];
#else /* weighted median */
switch(cnt)
{
case 0 :
Median[node] = -1;
hasfixed = 1;
break;
case 1 :
Median[node] = SCALE*Medadj[0];
break;
case 2 :
Median[node] = (Medadj[0]+Medadj[1])*(SCALE/2);
break;
default :
qsort((char*)Medadj,cnt,sizeof(Medadj[0]),(qsortcmpfun)intcmp);
if(cnt%2)
{
Median[node] = SCALE*Medadj[cnt/2];
break;
}
register int r = cnt/2;
register int l = r-1;
// compute left/right spans
int lspan = 0, rspan = 0;
int curinv = Medadj[l];
int nextinv;
for(k = l-1; k >= 0; --k)
{
nextinv = Medadj[k];
lspan += curinv-nextinv;
curinv = nextinv;
}
curinv = Medadj[r];
for(k = r+1; k < cnt; ++k)
{
nextinv = Medadj[k];
rspan += nextinv-curinv;
curinv = nextinv;
}
if(lspan == rspan)
Median[node] = (Medadj[l]+Medadj[r])*(SCALE/2);
else
{
int w = Medadj[l]*rspan+ Medadj[r]*lspan;
Median[node] = (w*SCALE)/(lspan+rspan);
}
break;
}
#endif /* RMEDIAN */
#endif /* BCENTER */
}
// sort by medians
if(mediansort(thislist,n_this,revert,hasfixed))
{
for(int i = 0; *thislist != -1; ++i)
Invert[*thislist++] = i;
if(thislev > 0)
Cross[thislev-1] = -1;
if(thislev < Maxlevel)
Cross[thislev] = -1;
}
}
for(lev = 0; lev < Maxlevel; ++lev)
if(Cross[lev] < 0)
Cross[lev] = rcross(ranks[lev],ranks[lev+1]);
}
// minimizing crossing using the median and transpose heuristics
static int mincross(int **ranks, int bestcross)
{
// copy current rank to List to start the process
copyrank(ranks,List);
inversion(List);
int curcross = cross_stat(List);
if(curcross < bestcross)
bestcross = curcross;
int counter = 1;
for (int iter = 0; iter < Maxiter; iter++)
{
median(List,iter%4 < 2 ? 1 : 0,iter%2 ? 1 : 0);
#ifndef NOTRANSPOSE
transpose(List,iter%4 < 2 ? 0 : 1);
#endif
curcross = 0;
for(int i = 0; i < Maxlevel; ++i)
curcross += Cross[i];
#ifdef BIGGRAPH
if(Verbose)
fprintf(stderr,"Iteration %d, curcross=%d bestcross=%d\n",
iter,curcross,bestcross);
#endif
if(curcross < bestcross)
{
float reduce = curcross <= 0 ? 0. :
((float)(bestcross-curcross))/curcross;
copyrank(List,ranks);
bestcross = curcross;
if(bestcross == 0 || reduce < Convergence)
break;
counter = 1;
}
else
{
if(counter > Minquit)
break;
++counter;
}
}
return bestcross;
}
// compute edge crossings
static int scross(int iw, int count, edge_t *adj, int dir)
{
int cross = 0;
if(dir > 0)
{
for(register edge_t *e = adj; e; e = e->next)
if(Invert[e->node] < iw)
cross += e->count*count;
}
else
{
for(register edge_t *e = adj; e; e = e->next)
if(Invert[e->node] > iw)
cross += e->count*count;
}
return cross;
}
// put stems in place where they'll be close to their adjacent nodes
static void fixstem(int **ranks, int *stem)
{
for(int dir = 1; dir >= -1; dir -= 2)
for(int l = dir > 0 ? 1 : Maxlevel-1; l >= 0 && l <= Maxlevel; l += dir)
{
edge_t **these = dir > 0 ? Inedges : Outedges;
edge_t **those = dir > 0 ? Outedges : Inedges;
int *rank = ranks[l];
// count the list size
int n_rank = 0;
int n_stem = 0;
for(int k = 0; rank[k] >= 0; ++k, ++n_rank)
{
int v = rank[k];
if(these[v] && !these[v]->next && !those[v])
stem[n_stem++] = v;
}
if(n_stem <= 0)
continue;
for(k = 0; k < n_stem; ++k)
{
int v = stem[k], iv = Invert[v];
int w = these[v]->node, iw = Invert[w];
int count = these[v]->count;
// find the left and right medians and bounds
int n = 0;
for(edge_t *e = those[w]; e; e = e->next)
if(e->node != v)
Medadj[n++] = Invert[e->node];
if(n == 0)
continue;
if(n > 2)
qsort((char*)Medadj,n,sizeof(Medadj[0]),(qsortcmpfun)intcmp);
int rmed = Medadj[n/2];
int lmed = n%2 ? rmed : Medadj[n/2-1];
// check each position
int best = iv, cross = 0;
int dist = iabs(2*iv - (rmed+lmed));
for(int d = -1; d <= 1; d += 2)
{
int c = 0, lm = lmed, rm = rmed;
for(int p = iv+d; p >= 0 && p < n_rank; p += d)
{
int x = rank[p];
if(p == lm)
lm -= d;
if(p == rm)
rm -= d;
int c0 = scross(iw,count,these[x],d);
int c1 = scross(iw,count,these[x],-d);
c += c1-c0;
if(c == cross)
{
int di = iabs(2*p -(lm+rm));
if(di < dist)
{
best = p;
dist = di;
}
}
else if(c < cross)
{
cross = c;
best = p;
dist = iabs(2*p - (lm+rm));
}
}
}
if(best == iv)
continue;
// move v to its best place
d = best < iv ? -1 : 1;
for(register int p = iv; p != best; p += d)
{
rank[p] = rank[p+d];
Invert[rank[p]] = p;
}
rank[p] = v;
Invert[v] = p;
}
}
}
// make node orderings for a connected component
static int buildrank(int **target, int **tryout, int component)
{
// build order using downward depth-first search
int start_target, start_tryout, bestcross, curcross;
initrank(target,component,1);
inversion(target);
start_target = bestcross = cross_stat(target);
start_tryout = 0;
if(bestcross > 0)
{
#ifndef NOTRANSPOSE
transpose(target,0);
bestcross = 0;
for(int i = 0; i < Maxlevel; ++i)
bestcross += Cross[i];
#endif
}
// rebuild order using upward depth-first search
if(bestcross > 0)
{
initrank(tryout,component,0);
inversion(tryout);
start_tryout = curcross = cross_stat(tryout);
if(curcross > 0)
{
#ifndef NOTRANSPOSE
transpose(tryout,0);
curcross = 0;
for(int i = 0; i < Maxlevel; ++i)
curcross += Cross[i];
#endif
}
if(curcross == 0)
{
copyrank(tryout,target);
bestcross = curcross;
}
}
// iterations to improve crossing
if(bestcross > 0)
bestcross = mincross(target,Maxint);
if(bestcross > 0 && (curcross = mincross(tryout,bestcross)) < bestcross)
copyrank(tryout,target);
// fix dangling stems
if(N_flat_edges <= 0)
{
inversion(target);
fixstem(target,Count);
}
return max(start_target,start_tryout);
}
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.