Actual source code: rcm.c
2: /* rcm.f -- translated by f2c (version 19931217).*/
4: #include <petscsys.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: /****************************************************************/
42: PetscErrorCode SPARSEPACKrcm(PetscInt *root, PetscInt *xadj, PetscInt *adjncy,
43: PetscInt *mask, PetscInt *perm, PetscInt *ccsize, PetscInt *deg)
44: {
45: /* System generated locals */
46: PetscInt i__1, i__2;
48: /* Local variables */
49: PetscInt node, fnbr, lnbr, i, j, k, l, lperm, jstop, jstrt;
51: PetscInt *, PetscInt *, PetscInt *, PetscInt *);
52: PetscInt lbegin, lvlend, nbr;
54: /* FIND THE DEGREES OF THE NODES IN THE */
55: /* COMPONENT SPECIFIED BY MASK AND ../../... */
56: /* ------------------------------------- */
60: /* Parameter adjustments */
61: --deg;
62: --perm;
63: --mask;
64: --adjncy;
65: --xadj;
68: SPARSEPACKdegree(root, &xadj[1], &adjncy[1], &mask[1], °[1], ccsize, &perm[1]);
69: mask[*root] = 0;
70: if (*ccsize <= 1) {
71: return(0);
72: }
73: lvlend = 0;
74: lnbr = 1;
75: /* LBEGIN AND LVLEND POINT TO THE BEGINNING AND */
76: /* THE END OF THE CURRENT LEVEL RESPECTIVELY. */
77: L100:
78: lbegin = lvlend + 1;
79: lvlend = lnbr;
80: i__1 = lvlend;
81: for (i = lbegin; i <= i__1; ++i) {
82: /* FOR EACH NODE IN CURRENT LEVEL ... */
83: node = perm[i];
84: jstrt = xadj[node];
85: jstop = xadj[node + 1] - 1;
87: /* FIND THE UNNUMBERED NEIGHBORS OF NODE. */
88: /* FNBR AND LNBR POINT TO THE FIRST AND LAST */
89: /* UNNUMBERED NEIGHBORS RESPECTIVELY OF THE CURRENT */
90: /* NODE IN PERM. */
91: fnbr = lnbr + 1;
92: i__2 = jstop;
93: for (j = jstrt; j <= i__2; ++j) {
94: nbr = adjncy[j];
95: if (!mask[nbr]) {
96: goto L200;
97: }
98: ++lnbr;
99: mask[nbr] = 0;
100: perm[lnbr] = nbr;
101: L200:
102: ;
103: }
104: if (fnbr >= lnbr) {
105: goto L600;
106: }
107: /* SORT THE NEIGHBORS OF NODE IN INCREASING */
108: /* ORDER BY DEGREE. LINEAR INSERTION IS USED.*/
109: k = fnbr;
110: L300:
111: l = k;
112: ++k;
113: nbr = perm[k];
114: L400:
115: if (l < fnbr) {
116: goto L500;
117: }
118: lperm = perm[l];
119: if (deg[lperm] <= deg[nbr]) {
120: goto L500;
121: }
122: perm[l + 1] = lperm;
123: --l;
124: goto L400;
125: L500:
126: perm[l + 1] = nbr;
127: if (k < lnbr) {
128: goto L300;
129: }
130: L600:
131: ;
132: }
133: if (lnbr > lvlend) {
134: goto L100;
135: }
136: /* WE NOW HAVE THE CUTHILL MCKEE ORDERING.*/
137: /* REVERSE IT BELOW ...*/
138: k = *ccsize / 2;
139: l = *ccsize;
140: i__1 = k;
141: for (i = 1; i <= i__1; ++i) {
142: lperm = perm[l];
143: perm[l] = perm[i];
144: perm[i] = lperm;
145: --l;
146: }
147: return(0);
148: }