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