Actual source code: qmdupd.c
2: /* qmdupd.f -- translated by f2c (version 19931217).*/
4: #include <petscsys.h>
6: /******************************************************************/
7: /*********** QMDUPD ..... QUOT MIN DEG UPDATE ************/
8: /******************************************************************/
9: /******************************************************************/
11: /* PURPOSE - THIS ROUTINE PERFORMS DEGREE UPDATE FOR A SET*/
12: /* OF NODES IN THE MINIMUM DEGREE ALGORITHM.*/
14: /* INPUT PARAMETERS -*/
15: /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE.*/
16: /* (NLIST, LIST) - THE LIST OF NODES WHOSE DEGREE HAS TO*/
17: /* BE UPDATED.*/
19: /* UPDATED PARAMETERS -*/
20: /* DEG - THE DEGREE VECTOR.*/
21: /* QSIZE - SIZE OF INDISTINGUISHABLE SUPERNODES.*/
22: /* QLINK - LINKED LIST FOR INDISTINGUISHABLE NODES.*/
23: /* MARKER - USED TO MARK THOSE NODES IN REACH/NBRHD SETS.*/
25: /* WORKING PARAMETERS -*/
26: /* RCHSET - THE REACHABLE SET.*/
27: /* NBRHD - THE NEIGHBORHOOD SET.*/
29: /* PROGRAM SUBROUTINES -*/
30: /* QMDMRG.*/
31: /******************************************************************/
34: PetscErrorCode SPARSEPACKqmdupd(PetscInt *xadj, PetscInt *adjncy, PetscInt *nlist,
35: PetscInt *list, PetscInt *deg, PetscInt *qsize, PetscInt *qlink, PetscInt *
36: marker, PetscInt *rchset, PetscInt *nbrhd)
37: {
38: /* System generated locals */
39: PetscInt i__1, i__2;
41: /* Local variables */
42: PetscInt inhd, irch, node, mark, j, inode, nabor, jstop, jstrt, il;
44: PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *),
45: SPARSEPACKqmdmrg(PetscInt*, PetscInt *, PetscInt *, PetscInt *, PetscInt *,
46: PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *);
47: PetscInt nhdsze, rchsze, deg0, deg1;
49: /* FIND ALL ELIMINATED SUPERNODES THAT ARE ADJACENT*/
50: /* TO SOME NODES IN THE GIVEN LIST. PUT THEM INTO.*/
51: /* (NHDSZE, NBRHD). DEG0 CONTAINS THE NUMBER OF*/
52: /* NODES IN THE LIST.*/
56: /* Parameter adjustments */
57: --nbrhd;
58: --rchset;
59: --marker;
60: --qlink;
61: --qsize;
62: --deg;
63: --list;
64: --adjncy;
65: --xadj;
67: if (*nlist <= 0) {
68: return(0);
69: }
70: deg0 = 0;
71: nhdsze = 0;
72: i__1 = *nlist;
73: for (il = 1; il <= i__1; ++il) {
74: node = list[il];
75: deg0 += qsize[node];
76: jstrt = xadj[node];
77: jstop = xadj[node + 1] - 1;
78: i__2 = jstop;
79: for (j = jstrt; j <= i__2; ++j) {
80: nabor = adjncy[j];
81: if (marker[nabor] != 0 || deg[nabor] >= 0) {
82: goto L100;
83: }
84: marker[nabor] = -1;
85: ++nhdsze;
86: nbrhd[nhdsze] = nabor;
87: L100:
88: ;
89: }
90: }
91: /* MERGE INDISTINGUISHABLE NODES IN THE LIST BY*/
92: /* CALLING THE SUBROUTINE QMDMRG.*/
93: if (nhdsze > 0) {
94: SPARSEPACKqmdmrg(&xadj[1], &adjncy[1], °[1], &qsize[1], &qlink[1], &marker[
95: 1], °0, &nhdsze, &nbrhd[1], &rchset[1], &nbrhd[nhdsze + 1]);
96: }
97: /* FIND THE NEW DEGREES OF THE NODES THAT HAVE NOT BEEN*/
98: /* MERGED.*/
99: i__1 = *nlist;
100: for (il = 1; il <= i__1; ++il) {
101: node = list[il];
102: mark = marker[node];
103: if (mark > 1 || mark < 0) {
104: goto L600;
105: }
106: marker[node] = 2;
107: SPARSEPACKqmdrch(&node, &xadj[1], &adjncy[1], °[1], &marker[1], &rchsze, &
108: rchset[1], &nhdsze, &nbrhd[1]);
109: deg1 = deg0;
110: if (rchsze <= 0) {
111: goto L400;
112: }
113: i__2 = rchsze;
114: for (irch = 1; irch <= i__2; ++irch) {
115: inode = rchset[irch];
116: deg1 += qsize[inode];
117: marker[inode] = 0;
118: }
119: L400:
120: deg[node] = deg1 - 1;
121: if (nhdsze <= 0) {
122: goto L600;
123: }
124: i__2 = nhdsze;
125: for (inhd = 1; inhd <= i__2; ++inhd) {
126: inode = nbrhd[inhd];
127: marker[inode] = 0;
128: }
129: L600:
130: ;
131: }
132: return(0);
133: }