Actual source code: rcm.c
1: /*$Id: rcm.c,v 1.16 2001/03/23 23:22:51 balay Exp $*/
2: /* rcm.f -- translated by f2c (version 19931217).*/
4: #include petsc.h
6: /*****************************************************************/
7: /********* RCM ..... REVERSE CUTHILL-MCKEE ORDERING *******/
8: /*****************************************************************/
9: /* PURPOSE - RCM NUMBERS A CONNECTED COMPONENT SPECIFIED BY */
10: /* MASK AND ../../.., USING THE RCM ALGORITHM. */
11: /* THE NUMBERING IS TO BE STARTED AT THE NODE ../../... */
12: /* */
13: /* INPUT PARAMETERS - */
14: /* ../../.. - IS THE NODE THAT DEFINES THE CONNECTED */
15: /* COMPONENT AND IT IS USED AS THE STARTING */
16: /* NODE FOR THE RCM ORDERING. */
17: /* (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR FOR */
18: /* THE GRAPH. */
19: /* */
20: /* UPDATED PARAMETERS - */
21: /* MASK - ONLY THOSE NODES WITH NONZERO INPUT MASK */
22: /* VALUES ARE CONSIDERED BY THE ROUTINE. THE */
23: /* NODES NUMBERED BY RCM WILL HAVE THEIR */
24: /* MASK VALUES SET TO ZERO. */
25: /* */
26: /* OUTPUT PARAMETERS - */
27: /* PERM - WILL CONTAIN THE RCM ORDERING. */
28: /* CCSIZE - IS THE SIZE OF THE CONNECTED COMPONENT */
29: /* THAT HAS BEEN NUMBERED BY RCM. */
30: /* */
31: /* WORKING PARAMETER - */
32: /* DEG - IS A TEMPORARY VECTOR USED TO HOLD THE DEGREE */
33: /* OF THE NODES IN THE SECTION GRAPH SPECIFIED */
34: /* BY MASK AND ../../... */
35: /* */
36: /* PROGRAM SUBROUTINES - */
37: /* DEGREE. */
38: /* */
39: /****************************************************************/
40: int SPARSEPACKrcm(int *root, int *xadj, int *adjncy,
41: int *mask, int *perm, int *ccsize, int *deg)
42: {
43: /* System generated locals */
44: int i__1, i__2;
46: /* Local variables */
47: int node, fnbr, lnbr, i, j, k, l, lperm, jstop, jstrt;
48: EXTERN int SPARSEPACKdegree(int *, int *, int *,
49: int *, int *, int *, int *);
50: int lbegin, lvlend, nbr;
52: /* FIND THE DEGREES OF THE NODES IN THE */
53: /* COMPONENT SPECIFIED BY MASK AND ../../... */
54: /* ------------------------------------- */
58: /* Parameter adjustments */
59: --deg;
60: --perm;
61: --mask;
62: --adjncy;
63: --xadj;
66: SPARSEPACKdegree(root, &xadj[1], &adjncy[1], &mask[1], °[1], ccsize, &perm[1]);
67: mask[*root] = 0;
68: if (*ccsize <= 1) {
69: return(0);
70: }
71: lvlend = 0;
72: lnbr = 1;
73: /* LBEGIN AND LVLEND POINT TO THE BEGINNING AND */
74: /* THE END OF THE CURRENT LEVEL RESPECTIVELY. */
75: L100:
76: lbegin = lvlend + 1;
77: lvlend = lnbr;
78: i__1 = lvlend;
79: for (i = lbegin; i <= i__1; ++i) {
80: /* FOR EACH NODE IN CURRENT LEVEL ... */
81: node = perm[i];
82: jstrt = xadj[node];
83: jstop = xadj[node + 1] - 1;
85: /* FIND THE UNNUMBERED NEIGHBORS OF NODE. */
86: /* FNBR AND LNBR POINT TO THE FIRST AND LAST */
87: /* UNNUMBERED NEIGHBORS RESPECTIVELY OF THE CURRENT */
88: /* NODE IN PERM. */
89: fnbr = lnbr + 1;
90: i__2 = jstop;
91: for (j = jstrt; j <= i__2; ++j) {
92: nbr = adjncy[j];
93: if (mask[nbr] == 0) {
94: goto L200;
95: }
96: ++lnbr;
97: mask[nbr] = 0;
98: perm[lnbr] = nbr;
99: L200:
100: ;
101: }
102: if (fnbr >= lnbr) {
103: goto L600;
104: }
105: /* SORT THE NEIGHBORS OF NODE IN INCREASING */
106: /* ORDER BY DEGREE. LINEAR INSERTION IS USED.*/
107: k = fnbr;
108: L300:
109: l = k;
110: ++k;
111: nbr = perm[k];
112: L400:
113: if (l < fnbr) {
114: goto L500;
115: }
116: lperm = perm[l];
117: if (deg[lperm] <= deg[nbr]) {
118: goto L500;
119: }
120: perm[l + 1] = lperm;
121: --l;
122: goto L400;
123: L500:
124: perm[l + 1] = nbr;
125: if (k < lnbr) {
126: goto L300;
127: }
128: L600:
129: ;
130: }
131: if (lnbr > lvlend) {
132: goto L100;
133: }
134: /* WE NOW HAVE THE CUTHILL MCKEE ORDERING.*/
135: /* REVERSE IT BELOW ...*/
136: k = *ccsize / 2;
137: l = *ccsize;
138: i__1 = k;
139: for (i = 1; i <= i__1; ++i) {
140: lperm = perm[l];
141: perm[l] = perm[i];
142: perm[i] = lperm;
143: --l;
144: }
145: return(0);
146: }