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