Actual source code: genqmd.c
1: /* genqmd.f -- translated by f2c (version 19931217).*/
3: #include petsc.h
5: /******************************************************************/
6: /*********** GENQMD ..... QUOT MIN DEGREE ORDERING **********/
7: /******************************************************************/
8: /* PURPOSE - THIS ROUTINE IMPLEMENTS THE MINIMUM DEGREE */
9: /* ALGORITHM. IT MAKES USE OF THE IMPLICIT REPRESENT- */
10: /* ATION OF THE ELIMINATION GRAPHS BY QUOTIENT GRAPHS, */
11: /* AND THE NOTION OF INDISTINGUISHABLE NODES. */
12: /* CAUTION - THE ADJACENCY VECTOR ADJNCY WILL BE */
13: /* DESTROYED. */
14: /* */
15: /* INPUT PARAMETERS - */
16: /* NEQNS - NUMBER OF EQUATIONS. */
17: /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE. */
18: /* */
19: /* OUTPUT PARAMETERS - */
20: /* PERM - THE MINIMUM DEGREE ORDERING. */
21: /* INVP - THE INVERSE OF PERM. */
22: /* */
23: /* WORKING PARAMETERS - */
24: /* DEG - THE DEGREE VECTOR. DEG(I) IS NEGATIVE MEANS */
25: /* NODE I HAS BEEN NUMBERED. */
26: /* MARKER - A MARKER VECTOR, WHERE MARKER(I) IS */
27: /* NEGATIVE MEANS NODE I HAS BEEN MERGED WITH */
28: /* ANOTHER NODE AND THUS CAN BE IGNORED. */
29: /* RCHSET - VECTOR USED FOR THE REACHABLE SET. */
30: /* NBRHD - VECTOR USED FOR THE NEIGHBORHOOD SET. */
31: /* QSIZE - VECTOR USED TO STORE THE SIZE OF */
32: /* INDISTINGUISHABLE SUPERNODES. */
33: /* QLINK - VECTOR TO STORE INDISTINGUISHABLE NODES, */
34: /* I, QLINK(I), QLINK(QLINK(I)) ... ARE THE */
35: /* MEMBERS OF THE SUPERNODE REPRESENTED BY I. */
36: /* */
37: /* PROGRAM SUBROUTINES - */
38: /* QMDRCH, QMDQT, QMDUPD. */
39: /* */
40: /******************************************************************/
41: /* */
42: /* */
45: PetscErrorCode SPARSEPACKgenqmd(PetscInt *neqns, PetscInt *xadj, PetscInt *adjncy,
46: PetscInt *perm, PetscInt *invp, PetscInt *deg, PetscInt *marker, PetscInt *
47: rchset, PetscInt *nbrhd, PetscInt *qsize, PetscInt *qlink, PetscInt *nofsub)
48: {
49: /* System generated locals */
50: PetscInt i__1;
52: /* Local variables */
53: PetscInt ndeg, irch, node, nump1, j, inode;
54: EXTERN PetscErrorCode SPARSEPACKqmdqt(PetscInt*, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *);
55: PetscInt ip, np, mindeg, search;
56: EXTERN PetscErrorCode SPARSEPACKqmdrch(PetscInt*, PetscInt *, PetscInt *,
57: PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *),
58: SPARSEPACKqmdupd(PetscInt*, PetscInt *, PetscInt *, PetscInt *, PetscInt *,
59: PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *);
60: PetscInt nhdsze, nxnode, rchsze, thresh, num;
62: /* INITIALIZE DEGREE VECTOR AND OTHER WORKING VARIABLES. */
65: /* Parameter adjustments */
66: --qlink;
67: --qsize;
68: --nbrhd;
69: --rchset;
70: --marker;
71: --deg;
72: --invp;
73: --perm;
74: --adjncy;
75: --xadj;
77: mindeg = *neqns;
78: *nofsub = 0;
79: i__1 = *neqns;
80: for (node = 1; node <= i__1; ++node) {
81: perm[node] = node;
82: invp[node] = node;
83: marker[node] = 0;
84: qsize[node] = 1;
85: qlink[node] = 0;
86: ndeg = xadj[node + 1] - xadj[node];
87: deg[node] = ndeg;
88: if (ndeg < mindeg) {
89: mindeg = ndeg;
90: }
91: }
92: num = 0;
93: /* PERFORM THRESHOLD SEARCH TO GET A NODE OF MIN DEGREE. */
94: /* VARIABLE SEARCH POINTS TO WHERE SEARCH SHOULD START. */
95: L200:
96: search = 1;
97: thresh = mindeg;
98: mindeg = *neqns;
99: L300:
100: nump1 = num + 1;
101: if (nump1 > search) {
102: search = nump1;
103: }
104: i__1 = *neqns;
105: for (j = search; j <= i__1; ++j) {
106: node = perm[j];
107: if (marker[node] < 0) {
108: goto L400;
109: }
110: ndeg = deg[node];
111: if (ndeg <= thresh) {
112: goto L500;
113: }
114: if (ndeg < mindeg) {
115: mindeg = ndeg;
116: }
117: L400:
118: ;
119: }
120: goto L200;
121: /* NODE HAS MINIMUM DEGREE. FIND ITS REACHABLE SETS BY */
122: /* CALLING QMDRCH. */
123: L500:
124: search = j;
125: *nofsub += deg[node];
126: marker[node] = 1;
127: SPARSEPACKqmdrch(&node, &xadj[1], &adjncy[1], °[1], &marker[1], &rchsze, &
128: rchset[1], &nhdsze, &nbrhd[1]);
129: /* ELIMINATE ALL NODES INDISTINGUISHABLE FROM NODE. */
130: /* THEY ARE GIVEN BY NODE, QLINK(NODE), .... */
131: nxnode = node;
132: L600:
133: ++num;
134: np = invp[nxnode];
135: ip = perm[num];
136: perm[np] = ip;
137: invp[ip] = np;
138: perm[num] = nxnode;
139: invp[nxnode] = num;
140: deg[nxnode] = -1;
141: nxnode = qlink[nxnode];
142: if (nxnode > 0) {
143: goto L600;
144: }
145: if (rchsze <= 0) {
146: goto L800;
147: }
148: /* UPDATE THE DEGREES OF THE NODES IN THE REACHABLE */
149: /* SET AND IDENTIFY INDISTINGUISHABLE NODES. */
150: SPARSEPACKqmdupd(&xadj[1], &adjncy[1], &rchsze, &rchset[1], °[1], &qsize[1], &
151: qlink[1], &marker[1], &rchset[rchsze + 1], &nbrhd[nhdsze + 1]);
152: /* RESET MARKER VALUE OF NODES IN REACH SET. */
153: /* UPDATE THRESHOLD VALUE FOR CYCLIC SEARCH. */
154: /* ALSO CALL QMDQT TO FORM NEW QUOTIENT GRAPH. */
155: marker[node] = 0;
156: i__1 = rchsze;
157: for (irch = 1; irch <= i__1; ++irch) {
158: inode = rchset[irch];
159: if (marker[inode] < 0) {
160: goto L700;
161: }
162: marker[inode] = 0;
163: ndeg = deg[inode];
164: if (ndeg < mindeg) {
165: mindeg = ndeg;
166: }
167: if (ndeg > thresh) {
168: goto L700;
169: }
170: mindeg = thresh;
171: thresh = ndeg;
172: search = invp[inode];
173: L700:
174: ;
175: }
176: if (nhdsze > 0) {
177: SPARSEPACKqmdqt(&node, &xadj[1], &adjncy[1], &marker[1], &rchsze, &rchset[1], &
178: nbrhd[1]);
179: }
180: L800:
181: if (num < *neqns) {
182: goto L300;
183: }
184: return(0);
185: }