Actual source code: genqmd.c
petsc-3.9.4 2018-09-11
2: /* genqmd.f -- translated by f2c (version 19931217).*/
4: #include <petscsys.h>
5: #include <petsc/private/matorderimpl.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: /* */
45: PetscErrorCode SPARSEPACKgenqmd(const PetscInt *neqns,const PetscInt *xadj,const PetscInt *adjncy,
46: PetscInt *perm, PetscInt *invp, PetscInt *deg, PetscInt *marker,
47: PetscInt *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: PetscInt ip, np, mindeg, search;
55: PetscInt nhdsze, nxnode, rchsze, thresh, num;
57: /* INITIALIZE DEGREE VECTOR AND OTHER WORKING VARIABLES. */
60: /* Parameter adjustments */
61: --qlink;
62: --qsize;
63: --nbrhd;
64: --rchset;
65: --marker;
66: --deg;
67: --invp;
68: --perm;
69: --adjncy;
70: --xadj;
72: mindeg = *neqns;
73: *nofsub = 0;
74: i__1 = *neqns;
75: for (node = 1; node <= i__1; ++node) {
76: perm[node] = node;
77: invp[node] = node;
78: marker[node] = 0;
79: qsize[node] = 1;
80: qlink[node] = 0;
81: ndeg = xadj[node + 1] - xadj[node];
82: deg[node] = ndeg;
83: if (ndeg < mindeg) mindeg = ndeg;
84: }
85: num = 0;
86: /* PERFORM THRESHOLD SEARCH TO GET A NODE OF MIN DEGREE. */
87: /* VARIABLE SEARCH POINTS TO WHERE SEARCH SHOULD START. */
88: L200:
89: search = 1;
90: thresh = mindeg;
91: mindeg = *neqns;
92: L300:
93: nump1 = num + 1;
94: if (nump1 > search) search = nump1;
95: i__1 = *neqns;
96: for (j = search; j <= i__1; ++j) {
97: node = perm[j];
98: if (marker[node] < 0) goto L400;
99: ndeg = deg[node];
100: if (ndeg <= thresh) goto L500;
101: if (ndeg < mindeg) mindeg = ndeg;
102: L400:
103: ;
104: }
105: goto L200;
106: /* NODE HAS MINIMUM DEGREE. FIND ITS REACHABLE SETS BY */
107: /* CALLING QMDRCH. */
108: L500:
109: search = j;
110: *nofsub += deg[node];
111: marker[node] = 1;
112: SPARSEPACKqmdrch(&node, &xadj[1], &adjncy[1], °[1], &marker[1], &rchsze, &
113: rchset[1], &nhdsze, &nbrhd[1]);
114: /* ELIMINATE ALL NODES INDISTINGUISHABLE FROM NODE. */
115: /* THEY ARE GIVEN BY NODE, QLINK(NODE), .... */
116: nxnode = node;
117: L600:
118: ++num;
119: np = invp[nxnode];
120: ip = perm[num];
121: perm[np] = ip;
122: invp[ip] = np;
123: perm[num] = nxnode;
124: invp[nxnode] = num;
125: deg[nxnode] = -1;
126: nxnode = qlink[nxnode];
127: if (nxnode > 0) goto L600;
128: if (rchsze <= 0) goto L800;
130: /* UPDATE THE DEGREES OF THE NODES IN THE REACHABLE */
131: /* SET AND IDENTIFY INDISTINGUISHABLE NODES. */
132: SPARSEPACKqmdupd(&xadj[1], &adjncy[1], &rchsze, &rchset[1], °[1], &qsize[1], &
133: qlink[1], &marker[1], &rchset[rchsze + 1], &nbrhd[nhdsze + 1]);
135: /* RESET MARKER VALUE OF NODES IN REACH SET. */
136: /* UPDATE THRESHOLD VALUE FOR CYCLIC SEARCH. */
137: /* ALSO CALL QMDQT TO FORM NEW QUOTIENT GRAPH. */
138: marker[node] = 0;
139: i__1 = rchsze;
140: for (irch = 1; irch <= i__1; ++irch) {
141: inode = rchset[irch];
142: if (marker[inode] < 0) goto L700;
144: marker[inode] = 0;
145: ndeg = deg[inode];
146: if (ndeg < mindeg) mindeg = ndeg;
147: if (ndeg > thresh) goto L700;
148: mindeg = thresh;
149: thresh = ndeg;
150: search = invp[inode];
151: L700:
152: ;
153: }
154: if (nhdsze > 0) {
155: SPARSEPACKqmdqt(&node, &xadj[1], &adjncy[1], &marker[1], &rchsze, &rchset[1], &
156: nbrhd[1]);
157: }
158: L800:
159: if (num < *neqns) goto L300;
160: return(0);
161: }