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