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