Actual source code: qmdmrg.c
1: #define PETSCMAT_DLL
3: /* qmdmrg.f -- translated by f2c (version 19931217).*/
5: #include petscsys.h
7: /******************************************************************/
8: /*********** QMDMRG ..... QUOT MIN DEG MERGE ************/
9: /******************************************************************/
10: /* PURPOSE - THIS ROUTINE MERGES INDISTINGUISHABLE NODES IN */
11: /* THE MINIMUM DEGREE ORDERING ALGORITHM. */
12: /* IT ALSO COMPUTES THE NEW DEGREES OF THESE */
13: /* NEW SUPERNODES. */
14: /* */
15: /* INPUT PARAMETERS - */
16: /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE. */
17: /* DEG0 - THE NUMBER OF NODES IN THE GIVEN SET. */
18: /* (NHDSZE, NBRHD) - THE SET OF ELIMINATED SUPERNODES */
19: /* ADJACENT TO SOME NODES IN THE SET. */
20: /* */
21: /* UPDATED PARAMETERS - */
22: /* DEG - THE DEGREE VECTOR. */
23: /* QSIZE - SIZE OF INDISTINGUISHABLE NODES. */
24: /* QLINK - LINKED LIST FOR INDISTINGUISHABLE NODES. */
25: /* MARKER - THE GIVEN SET IS GIVEN BY THOSE NODES WITH */
26: /* MARKER VALUE SET TO 1. THOSE NODES WITH DEGREE */
27: /* UPDATED WILL HAVE MARKER VALUE SET TO 2. */
28: /* */
29: /* WORKING PARAMETERS - */
30: /* RCHSET - THE REACHABLE SET. */
31: /* OVRLP - TEMP VECTOR TO STORE THE INTERSECTION OF TWO */
32: /* REACHABLE SETS. */
33: /* */
34: /*****************************************************************/
37: PetscErrorCode SPARSEPACKqmdmrg(PetscInt *xadj, PetscInt *adjncy, PetscInt *deg,
38: PetscInt *qsize, PetscInt *qlink, PetscInt *marker, PetscInt *deg0,
39: PetscInt *nhdsze, PetscInt *nbrhd, PetscInt *rchset, PetscInt *ovrlp)
40: {
41: /* System generated locals */
42: PetscInt i__1, i__2, i__3;
44: /* Local variables */
45: PetscInt head, inhd, irch, node, mark, ilink, root, j, lnode, nabor,
46: jstop, jstrt, rchsze, mrgsze, novrlp, iov, deg1;
49: /* Parameter adjustments */
50: --ovrlp;
51: --rchset;
52: --nbrhd;
53: --marker;
54: --qlink;
55: --qsize;
56: --deg;
57: --adjncy;
58: --xadj;
60: if (*nhdsze <= 0) {
61: return(0);
62: }
63: i__1 = *nhdsze;
64: for (inhd = 1; inhd <= i__1; ++inhd) {
65: root = nbrhd[inhd];
66: marker[root] = 0;
67: }
68: /* LOOP THROUGH EACH ELIMINATED SUPERNODE IN THE SET */
69: /* (NHDSZE, NBRHD). */
70: i__1 = *nhdsze;
71: for (inhd = 1; inhd <= i__1; ++inhd) {
72: root = nbrhd[inhd];
73: marker[root] = -1;
74: rchsze = 0;
75: novrlp = 0;
76: deg1 = 0;
77: L200:
78: jstrt = xadj[root];
79: jstop = xadj[root + 1] - 1;
80: /* DETERMINE THE REACHABLE SET AND ITS PETSCINTERSECT- */
81: /* ION WITH THE INPUT REACHABLE SET. */
82: i__2 = jstop;
83: for (j = jstrt; j <= i__2; ++j) {
84: nabor = adjncy[j];
85: root = -nabor;
86: if (nabor < 0) {
87: goto L200;
88: } else if (!nabor) {
89: goto L700;
90: } else {
91: goto L300;
92: }
93: L300:
94: mark = marker[nabor];
95: if (mark < 0) {
96: goto L600;
97: } else if (!mark) {
98: goto L400;
99: } else {
100: goto L500;
101: }
102: L400:
103: ++rchsze;
104: rchset[rchsze] = nabor;
105: deg1 += qsize[nabor];
106: marker[nabor] = 1;
107: goto L600;
108: L500:
109: if (mark > 1) {
110: goto L600;
111: }
112: ++novrlp;
113: ovrlp[novrlp] = nabor;
114: marker[nabor] = 2;
115: L600:
116: ;
117: }
118: /* FROM THE OVERLAPPED SET, DETERMINE THE NODES */
119: /* THAT CAN BE MERGED TOGETHER. */
120: L700:
121: head = 0;
122: mrgsze = 0;
123: i__2 = novrlp;
124: for (iov = 1; iov <= i__2; ++iov) {
125: node = ovrlp[iov];
126: jstrt = xadj[node];
127: jstop = xadj[node + 1] - 1;
128: i__3 = jstop;
129: for (j = jstrt; j <= i__3; ++j) {
130: nabor = adjncy[j];
131: if (marker[nabor] != 0) {
132: goto L800;
133: }
134: marker[node] = 1;
135: goto L1100;
136: L800:
137: ;
138: }
139: /* NODE BELONGS TO THE NEW MERGED SUPERNODE. */
140: /* UPDATE THE VECTORS QLINK AND QSIZE. */
141: mrgsze += qsize[node];
142: marker[node] = -1;
143: lnode = node;
144: L900:
145: ilink = qlink[lnode];
146: if (ilink <= 0) {
147: goto L1000;
148: }
149: lnode = ilink;
150: goto L900;
151: L1000:
152: qlink[lnode] = head;
153: head = node;
154: L1100:
155: ;
156: }
157: if (head <= 0) {
158: goto L1200;
159: }
160: qsize[head] = mrgsze;
161: deg[head] = *deg0 + deg1 - 1;
162: marker[head] = 2;
163: /* RESET MARKER VALUES. */
164: L1200:
165: root = nbrhd[inhd];
166: marker[root] = 0;
167: if (rchsze <= 0) {
168: goto L1400;
169: }
170: i__2 = rchsze;
171: for (irch = 1; irch <= i__2; ++irch) {
172: node = rchset[irch];
173: marker[node] = 0;
174: }
175: L1400:
176: ;
177: }
178: return(0);
179: }