1 /* qmdupd.f -- translated by f2c (version 19931217).*/
2
3 #include <petscsys.h>
4 #include <petsc/private/matorderimpl.h>
5
6 /******************************************************************/
7 /*********** QMDUPD ..... QUOT MIN DEG UPDATE ************/
8 /******************************************************************/
9 /******************************************************************/
10
11 /* PURPOSE - THIS ROUTINE PERFORMS DEGREE UPDATE FOR A SET*/
12 /* OF NODES IN THE MINIMUM DEGREE ALGORITHM.*/
13
14 /* INPUT PARAMETERS -*/
15 /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE.*/
16 /* (NLIST, LIST) - THE LIST OF NODES WHOSE DEGREE HAS TO*/
17 /* BE UPDATED.*/
18
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.*/
24
25 /* WORKING PARAMETERS -*/
26 /* RCHSET - THE REACHABLE SET.*/
27 /* NBRHD - THE NEIGHBORHOOD SET.*/
28
29 /* PROGRAM SUBROUTINES -*/
30 /* QMDMRG.*/
31 /******************************************************************/
SPARSEPACKqmdupd(const PetscInt * xadj,const PetscInt * adjncy,const PetscInt * nlist,const PetscInt * list,PetscInt * deg,PetscInt * qsize,PetscInt * qlink,PetscInt * marker,PetscInt * rchset,PetscInt * nbrhd)32 PetscErrorCode SPARSEPACKqmdupd(const PetscInt *xadj, const PetscInt *adjncy, const PetscInt *nlist, const PetscInt *list, PetscInt *deg, PetscInt *qsize, PetscInt *qlink, PetscInt *marker, PetscInt *rchset, PetscInt *nbrhd)
33 {
34 /* System generated locals */
35 PetscInt i__1, i__2;
36
37 /* Local variables */
38 PetscInt inhd, irch, node, mark, j, inode, nabor, jstop, jstrt, il;
39 PetscInt nhdsze, rchsze, deg0, deg1;
40
41 /* FIND ALL ELIMINATED SUPERNODES THAT ARE ADJACENT*/
42 /* TO SOME NODES IN THE GIVEN LIST. PUT THEM INTO.*/
43 /* (NHDSZE, NBRHD). DEG0 CONTAINS THE NUMBER OF*/
44 /* NODES IN THE LIST.*/
45
46 PetscFunctionBegin;
47 /* Parameter adjustments */
48 --nbrhd;
49 --rchset;
50 --marker;
51 --qlink;
52 --qsize;
53 --deg;
54 --list;
55 --adjncy;
56 --xadj;
57
58 if (*nlist <= 0) PetscFunctionReturn(PETSC_SUCCESS);
59 deg0 = 0;
60 nhdsze = 0;
61 i__1 = *nlist;
62 for (il = 1; il <= i__1; ++il) {
63 node = list[il];
64 deg0 += qsize[node];
65 jstrt = xadj[node];
66 jstop = xadj[node + 1] - 1;
67 i__2 = jstop;
68 for (j = jstrt; j <= i__2; ++j) {
69 nabor = adjncy[j];
70 if (marker[nabor] != 0 || deg[nabor] >= 0) goto L100;
71 marker[nabor] = -1;
72 ++nhdsze;
73 nbrhd[nhdsze] = nabor;
74 L100:;
75 }
76 }
77 /* MERGE INDISTINGUISHABLE NODES IN THE LIST BY*/
78 /* CALLING THE SUBROUTINE QMDMRG.*/
79 if (nhdsze > 0) PetscCall(SPARSEPACKqmdmrg(&xadj[1], &adjncy[1], °[1], &qsize[1], &qlink[1], &marker[1], °0, &nhdsze, &nbrhd[1], &rchset[1], &nbrhd[nhdsze + 1]));
80 /* FIND THE NEW DEGREES OF THE NODES THAT HAVE NOT BEEN*/
81 /* MERGED.*/
82 i__1 = *nlist;
83 for (il = 1; il <= i__1; ++il) {
84 node = list[il];
85 mark = marker[node];
86 if (mark > 1 || mark < 0) goto L600;
87 marker[node] = 2;
88 PetscCall(SPARSEPACKqmdrch(&node, &xadj[1], &adjncy[1], °[1], &marker[1], &rchsze, &rchset[1], &nhdsze, &nbrhd[1]));
89 deg1 = deg0;
90 if (rchsze <= 0) goto L400;
91 i__2 = rchsze;
92 for (irch = 1; irch <= i__2; ++irch) {
93 inode = rchset[irch];
94 deg1 += qsize[inode];
95 marker[inode] = 0;
96 }
97 L400:
98 deg[node] = deg1 - 1;
99 if (nhdsze <= 0) goto L600;
100 i__2 = nhdsze;
101 for (inhd = 1; inhd <= i__2; ++inhd) {
102 inode = nbrhd[inhd];
103 marker[inode] = 0;
104 }
105 L600:;
106 }
107 PetscFunctionReturn(PETSC_SUCCESS);
108 }
109