diff options
Diffstat (limited to 'benchmarks/CUDA/DG/3rdParty/ParMetis-3.1/METISLib/mmd.c')
| -rw-r--r-- | benchmarks/CUDA/DG/3rdParty/ParMetis-3.1/METISLib/mmd.c | 593 |
1 files changed, 593 insertions, 0 deletions
diff --git a/benchmarks/CUDA/DG/3rdParty/ParMetis-3.1/METISLib/mmd.c b/benchmarks/CUDA/DG/3rdParty/ParMetis-3.1/METISLib/mmd.c new file mode 100644 index 0000000..1b43618 --- /dev/null +++ b/benchmarks/CUDA/DG/3rdParty/ParMetis-3.1/METISLib/mmd.c @@ -0,0 +1,593 @@ +/* + * mmd.c + * + * ************************************************************** + * The following C function was developed from a FORTRAN subroutine + * in SPARSPAK written by Eleanor Chu, Alan George, Joseph Liu + * and Esmond Ng. + * + * The FORTRAN-to-C transformation and modifications such as dynamic + * memory allocation and deallocation were performed by Chunguang + * Sun. + * ************************************************************** + * + * Taken from SMMS, George 12/13/94 + * + * The meaning of invperm, and perm vectors is different from that + * in genqmd_ of SparsPak + * + * $Id: mmd.c,v 1.1 2003/07/16 15:55:11 karypis Exp $ + */ + +#include <metis.h> + + +/************************************************************************* +* genmmd -- multiple minimum external degree +* purpose -- this routine implements the minimum degree +* algorithm. it makes use of the implicit representation +* of elimination graphs by quotient graphs, and the notion +* of indistinguishable nodes. It also implements the modifications +* by multiple elimination and minimum external degree. +* Caution -- the adjacency vector adjncy will be destroyed. +* Input parameters -- +* neqns -- number of equations. +* (xadj, adjncy) -- the adjacency structure. +* delta -- tolerance value for multiple elimination. +* maxint -- maximum machine representable (short) integer +* (any smaller estimate will do) for marking nodes. +* Output parameters -- +* perm -- the minimum degree ordering. +* invp -- the inverse of perm. +* *ncsub -- an upper bound on the number of nonzero subscripts +* for the compressed storage scheme. +* Working parameters -- +* head -- vector for head of degree lists. +* invp -- used temporarily for degree forward link. +* perm -- used temporarily for degree backward link. +* qsize -- vector for size of supernodes. +* list -- vector for temporary linked lists. +* marker -- a temporary marker vector. +* Subroutines used -- mmdelm, mmdint, mmdnum, mmdupd. +**************************************************************************/ +void genmmd(int neqns, idxtype *xadj, idxtype *adjncy, idxtype *invp, idxtype *perm, + int delta, idxtype *head, idxtype *qsize, idxtype *list, idxtype *marker, + int maxint, int *ncsub) +{ + int ehead, i, mdeg, mdlmt, mdeg_node, nextmd, num, tag; + + if (neqns <= 0) + return; + + /* Adjust from C to Fortran */ + xadj--; adjncy--; invp--; perm--; head--; qsize--; list--; marker--; + + /* initialization for the minimum degree algorithm. */ + *ncsub = 0; + mmdint(neqns, xadj, adjncy, head, invp, perm, qsize, list, marker); + + /* 'num' counts the number of ordered nodes plus 1. */ + num = 1; + + /* eliminate all isolated nodes. */ + nextmd = head[1]; + while (nextmd > 0) { + mdeg_node = nextmd; + nextmd = invp[mdeg_node]; + marker[mdeg_node] = maxint; + invp[mdeg_node] = -num; + num = num + 1; + } + + /* search for node of the minimum degree. 'mdeg' is the current */ + /* minimum degree; 'tag' is used to facilitate marking nodes. */ + if (num > neqns) + goto n1000; + tag = 1; + head[1] = 0; + mdeg = 2; + + /* infinite loop here ! */ + while (1) { + while (head[mdeg] <= 0) + mdeg++; + + /* use value of 'delta' to set up 'mdlmt', which governs */ + /* when a degree update is to be performed. */ + mdlmt = mdeg + delta; + ehead = 0; + +n500: + mdeg_node = head[mdeg]; + while (mdeg_node <= 0) { + mdeg++; + + if (mdeg > mdlmt) + goto n900; + mdeg_node = head[mdeg]; + }; + + /* remove 'mdeg_node' from the degree structure. */ + nextmd = invp[mdeg_node]; + head[mdeg] = nextmd; + if (nextmd > 0) + perm[nextmd] = -mdeg; + invp[mdeg_node] = -num; + *ncsub += mdeg + qsize[mdeg_node] - 2; + if ((num+qsize[mdeg_node]) > neqns) + goto n1000; + + /* eliminate 'mdeg_node' and perform quotient graph */ + /* transformation. reset 'tag' value if necessary. */ + tag++; + if (tag >= maxint) { + tag = 1; + for (i = 1; i <= neqns; i++) + if (marker[i] < maxint) + marker[i] = 0; + }; + + mmdelm(mdeg_node, xadj, adjncy, head, invp, perm, qsize, list, marker, maxint, tag); + + num += qsize[mdeg_node]; + list[mdeg_node] = ehead; + ehead = mdeg_node; + if (delta >= 0) + goto n500; + + n900: + /* update degrees of the nodes involved in the */ + /* minimum degree nodes elimination. */ + if (num > neqns) + goto n1000; + mmdupd( ehead, neqns, xadj, adjncy, delta, &mdeg, head, invp, perm, qsize, list, marker, maxint, &tag); + }; /* end of -- while ( 1 ) -- */ + +n1000: + mmdnum( neqns, perm, invp, qsize ); + + /* Adjust from Fortran back to C*/ + xadj++; adjncy++; invp++; perm++; head++; qsize++; list++; marker++; +} + + +/************************************************************************** +* mmdelm ...... multiple minimum degree elimination +* Purpose -- This routine eliminates the node mdeg_node of minimum degree +* from the adjacency structure, which is stored in the quotient +* graph format. It also transforms the quotient graph representation +* of the elimination graph. +* Input parameters -- +* mdeg_node -- node of minimum degree. +* maxint -- estimate of maximum representable (short) integer. +* tag -- tag value. +* Updated parameters -- +* (xadj, adjncy) -- updated adjacency structure. +* (head, forward, backward) -- degree doubly linked structure. +* qsize -- size of supernode. +* marker -- marker vector. +* list -- temporary linked list of eliminated nabors. +***************************************************************************/ +void mmdelm(int mdeg_node, idxtype *xadj, idxtype *adjncy, idxtype *head, idxtype *forward, + idxtype *backward, idxtype *qsize, idxtype *list, idxtype *marker, int maxint,int tag) +{ + int element, i, istop, istart, j, + jstop, jstart, link, + nabor, node, npv, nqnbrs, nxnode, + pvnode, rlmt, rloc, rnode, xqnbr; + + /* find the reachable set of 'mdeg_node' and */ + /* place it in the data structure. */ + marker[mdeg_node] = tag; + istart = xadj[mdeg_node]; + istop = xadj[mdeg_node+1] - 1; + + /* 'element' points to the beginning of the list of */ + /* eliminated nabors of 'mdeg_node', and 'rloc' gives the */ + /* storage location for the next reachable node. */ + element = 0; + rloc = istart; + rlmt = istop; + for ( i = istart; i <= istop; i++ ) { + nabor = adjncy[i]; + if ( nabor == 0 ) break; + if ( marker[nabor] < tag ) { + marker[nabor] = tag; + if ( forward[nabor] < 0 ) { + list[nabor] = element; + element = nabor; + } else { + adjncy[rloc] = nabor; + rloc++; + }; + }; /* end of -- if -- */ + }; /* end of -- for -- */ + + /* merge with reachable nodes from generalized elements. */ + while ( element > 0 ) { + adjncy[rlmt] = -element; + link = element; + +n400: + jstart = xadj[link]; + jstop = xadj[link+1] - 1; + for ( j = jstart; j <= jstop; j++ ) { + node = adjncy[j]; + link = -node; + if ( node < 0 ) goto n400; + if ( node == 0 ) break; + if ((marker[node]<tag)&&(forward[node]>=0)) { + marker[node] = tag; + /*use storage from eliminated nodes if necessary.*/ + while ( rloc >= rlmt ) { + link = -adjncy[rlmt]; + rloc = xadj[link]; + rlmt = xadj[link+1] - 1; + }; + adjncy[rloc] = node; + rloc++; + }; + }; /* end of -- for ( j = jstart; -- */ + element = list[element]; + }; /* end of -- while ( element > 0 ) -- */ + if ( rloc <= rlmt ) adjncy[rloc] = 0; + /* for each node in the reachable set, do the following. */ + link = mdeg_node; + +n1100: + istart = xadj[link]; + istop = xadj[link+1] - 1; + for ( i = istart; i <= istop; i++ ) { + rnode = adjncy[i]; + link = -rnode; + if ( rnode < 0 ) goto n1100; + if ( rnode == 0 ) return; + + /* 'rnode' is in the degree list structure. */ + pvnode = backward[rnode]; + if (( pvnode != 0 ) && ( pvnode != (-maxint) )) { + /* then remove 'rnode' from the structure. */ + nxnode = forward[rnode]; + if ( nxnode > 0 ) backward[nxnode] = pvnode; + if ( pvnode > 0 ) forward[pvnode] = nxnode; + npv = -pvnode; + if ( pvnode < 0 ) head[npv] = nxnode; + }; + + /* purge inactive quotient nabors of 'rnode'. */ + jstart = xadj[rnode]; + jstop = xadj[rnode+1] - 1; + xqnbr = jstart; + for ( j = jstart; j <= jstop; j++ ) { + nabor = adjncy[j]; + if ( nabor == 0 ) break; + if ( marker[nabor] < tag ) { + adjncy[xqnbr] = nabor; + xqnbr++; + }; + }; + + /* no active nabor after the purging. */ + nqnbrs = xqnbr - jstart; + if ( nqnbrs <= 0 ) { + /* merge 'rnode' with 'mdeg_node'. */ + qsize[mdeg_node] += qsize[rnode]; + qsize[rnode] = 0; + marker[rnode] = maxint; + forward[rnode] = -mdeg_node; + backward[rnode] = -maxint; + } else { + /* flag 'rnode' for degree update, and */ + /* add 'mdeg_node' as a nabor of 'rnode'. */ + forward[rnode] = nqnbrs + 1; + backward[rnode] = 0; + adjncy[xqnbr] = mdeg_node; + xqnbr++; + if ( xqnbr <= jstop ) adjncy[xqnbr] = 0; + }; + }; /* end of -- for ( i = istart; -- */ + return; + } + +/*************************************************************************** +* mmdint ---- mult minimum degree initialization +* purpose -- this routine performs initialization for the +* multiple elimination version of the minimum degree algorithm. +* input parameters -- +* neqns -- number of equations. +* (xadj, adjncy) -- adjacency structure. +* output parameters -- +* (head, dfrow, backward) -- degree doubly linked structure. +* qsize -- size of supernode ( initialized to one). +* list -- linked list. +* marker -- marker vector. +****************************************************************************/ +int mmdint(int neqns, idxtype *xadj, idxtype *adjncy, idxtype *head, idxtype *forward, + idxtype *backward, idxtype *qsize, idxtype *list, idxtype *marker) +{ + int fnode, ndeg, node; + + for ( node = 1; node <= neqns; node++ ) { + head[node] = 0; + qsize[node] = 1; + marker[node] = 0; + list[node] = 0; + }; + + /* initialize the degree doubly linked lists. */ + for ( node = 1; node <= neqns; node++ ) { + ndeg = xadj[node+1] - xadj[node]/* + 1*/; /* george */ + if (ndeg == 0) + ndeg = 1; + fnode = head[ndeg]; + forward[node] = fnode; + head[ndeg] = node; + if ( fnode > 0 ) backward[fnode] = node; + backward[node] = -ndeg; + }; + return 0; +} + +/**************************************************************************** +* mmdnum --- multi minimum degree numbering +* purpose -- this routine performs the final step in producing +* the permutation and inverse permutation vectors in the +* multiple elimination version of the minimum degree +* ordering algorithm. +* input parameters -- +* neqns -- number of equations. +* qsize -- size of supernodes at elimination. +* updated parameters -- +* invp -- inverse permutation vector. on input, +* if qsize[node] = 0, then node has been merged +* into the node -invp[node]; otherwise, +* -invp[node] is its inverse labelling. +* output parameters -- +* perm -- the permutation vector. +****************************************************************************/ +void mmdnum(int neqns, idxtype *perm, idxtype *invp, idxtype *qsize) +{ + int father, nextf, node, nqsize, num, root; + + for ( node = 1; node <= neqns; node++ ) { + nqsize = qsize[node]; + if ( nqsize <= 0 ) perm[node] = invp[node]; + if ( nqsize > 0 ) perm[node] = -invp[node]; + }; + + /* for each node which has been merged, do the following. */ + for ( node = 1; node <= neqns; node++ ) { + if ( perm[node] <= 0 ) { + + /* trace the merged tree until one which has not */ + /* been merged, call it root. */ + father = node; + while ( perm[father] <= 0 ) + father = - perm[father]; + + /* number node after root. */ + root = father; + num = perm[root] + 1; + invp[node] = -num; + perm[root] = num; + + /* shorten the merged tree. */ + father = node; + nextf = - perm[father]; + while ( nextf > 0 ) { + perm[father] = -root; + father = nextf; + nextf = -perm[father]; + }; + }; /* end of -- if ( perm[node] <= 0 ) -- */ + }; /* end of -- for ( node = 1; -- */ + + /* ready to compute perm. */ + for ( node = 1; node <= neqns; node++ ) { + num = -invp[node]; + invp[node] = num; + perm[num] = node; + }; + return; +} + +/**************************************************************************** +* mmdupd ---- multiple minimum degree update +* purpose -- this routine updates the degrees of nodes after a +* multiple elimination step. +* input parameters -- +* ehead -- the beginning of the list of eliminated nodes +* (i.e., newly formed elements). +* neqns -- number of equations. +* (xadj, adjncy) -- adjacency structure. +* delta -- tolerance value for multiple elimination. +* maxint -- maximum machine representable (short) integer. +* updated parameters -- +* mdeg -- new minimum degree after degree update. +* (head, forward, backward) -- degree doubly linked structure. +* qsize -- size of supernode. +* list -- marker vector for degree update. +* *tag -- tag value. +****************************************************************************/ +void mmdupd(int ehead, int neqns, idxtype *xadj, idxtype *adjncy, int delta, int *mdeg, + idxtype *head, idxtype *forward, idxtype *backward, idxtype *qsize, idxtype *list, + idxtype *marker, int maxint,int *tag) +{ + int deg, deg0, element, enode, fnode, i, iq2, istop, + istart, j, jstop, jstart, link, mdeg0, mtag, nabor, + node, q2head, qxhead; + + mdeg0 = *mdeg + delta; + element = ehead; + +n100: + if ( element <= 0 ) return; + + /* for each of the newly formed element, do the following. */ + /* reset tag value if necessary. */ + mtag = *tag + mdeg0; + if ( mtag >= maxint ) { + *tag = 1; + for ( i = 1; i <= neqns; i++ ) + if ( marker[i] < maxint ) marker[i] = 0; + mtag = *tag + mdeg0; + }; + + /* create two linked lists from nodes associated with 'element': */ + /* one with two nabors (q2head) in the adjacency structure, and the*/ + /* other with more than two nabors (qxhead). also compute 'deg0',*/ + /* number of nodes in this element. */ + q2head = 0; + qxhead = 0; + deg0 = 0; + link =element; + +n400: + istart = xadj[link]; + istop = xadj[link+1] - 1; + for ( i = istart; i <= istop; i++ ) { + enode = adjncy[i]; + link = -enode; + if ( enode < 0 ) goto n400; + if ( enode == 0 ) break; + if ( qsize[enode] != 0 ) { + deg0 += qsize[enode]; + marker[enode] = mtag; + + /*'enode' requires a degree update*/ + if ( backward[enode] == 0 ) { + /* place either in qxhead or q2head list. */ + if ( forward[enode] != 2 ) { + list[enode] = qxhead; + qxhead = enode; + } else { + list[enode] = q2head; + q2head = enode; + }; + }; + }; /* enf of -- if ( qsize[enode] != 0 ) -- */ + }; /* end of -- for ( i = istart; -- */ + + /* for each node in q2 list, do the following. */ + enode = q2head; + iq2 = 1; + +n900: + if ( enode <= 0 ) goto n1500; + if ( backward[enode] != 0 ) goto n2200; + (*tag)++; + deg = deg0; + + /* identify the other adjacent element nabor. */ + istart = xadj[enode]; + nabor = adjncy[istart]; + if ( nabor == element ) nabor = adjncy[istart+1]; + link = nabor; + if ( forward[nabor] >= 0 ) { + /* nabor is uneliminated, increase degree count. */ + deg += qsize[nabor]; + goto n2100; + }; + + /* the nabor is eliminated. for each node in the 2nd element */ + /* do the following. */ +n1000: + istart = xadj[link]; + istop = xadj[link+1] - 1; + for ( i = istart; i <= istop; i++ ) { + node = adjncy[i]; + link = -node; + if ( node != enode ) { + if ( node < 0 ) goto n1000; + if ( node == 0 ) goto n2100; + if ( qsize[node] != 0 ) { + if ( marker[node] < *tag ) { + /* 'node' is not yet considered. */ + marker[node] = *tag; + deg += qsize[node]; + } else { + if ( backward[node] == 0 ) { + if ( forward[node] == 2 ) { + /* 'node' is indistinguishable from 'enode'.*/ + /* merge them into a new supernode. */ + qsize[enode] += qsize[node]; + qsize[node] = 0; + marker[node] = maxint; + forward[node] = -enode; + backward[node] = -maxint; + } else { + /* 'node' is outmacthed by 'enode' */ + if (backward[node]==0) backward[node] = -maxint; + }; + }; /* end of -- if ( backward[node] == 0 ) -- */ + }; /* end of -- if ( marker[node] < *tag ) -- */ + }; /* end of -- if ( qsize[node] != 0 ) -- */ + }; /* end of -- if ( node != enode ) -- */ + }; /* end of -- for ( i = istart; -- */ + goto n2100; + +n1500: + /* for each 'enode' in the 'qx' list, do the following. */ + enode = qxhead; + iq2 = 0; + +n1600: if ( enode <= 0 ) goto n2300; + if ( backward[enode] != 0 ) goto n2200; + (*tag)++; + deg = deg0; + + /*for each unmarked nabor of 'enode', do the following.*/ + istart = xadj[enode]; + istop = xadj[enode+1] - 1; + for ( i = istart; i <= istop; i++ ) { + nabor = adjncy[i]; + if ( nabor == 0 ) break; + if ( marker[nabor] < *tag ) { + marker[nabor] = *tag; + link = nabor; + if ( forward[nabor] >= 0 ) + /*if uneliminated, include it in deg count.*/ + deg += qsize[nabor]; + else { +n1700: + /* if eliminated, include unmarked nodes in this*/ + /* element into the degree count. */ + jstart = xadj[link]; + jstop = xadj[link+1] - 1; + for ( j = jstart; j <= jstop; j++ ) { + node = adjncy[j]; + link = -node; + if ( node < 0 ) goto n1700; + if ( node == 0 ) break; + if ( marker[node] < *tag ) { + marker[node] = *tag; + deg += qsize[node]; + }; + }; /* end of -- for ( j = jstart; -- */ + }; /* end of -- if ( forward[nabor] >= 0 ) -- */ + }; /* end of -- if ( marker[nabor] < *tag ) -- */ + }; /* end of -- for ( i = istart; -- */ + +n2100: + /* update external degree of 'enode' in degree structure, */ + /* and '*mdeg' if necessary. */ + deg = deg - qsize[enode] + 1; + fnode = head[deg]; + forward[enode] = fnode; + backward[enode] = -deg; + if ( fnode > 0 ) backward[fnode] = enode; + head[deg] = enode; + if ( deg < *mdeg ) *mdeg = deg; + +n2200: + /* get next enode in current element. */ + enode = list[enode]; + if ( iq2 == 1 ) goto n900; + goto n1600; + +n2300: + /* get next element in the list. */ + *tag = mtag; + element = list[element]; + goto n100; + } |