1 /* qmdqt.f -- translated by f2c (version 19931217).*/ 2 3 #include <petscsys.h> 4 #include <petsc/private/matorderimpl.h> 5 6 /***************************************************************/ 7 /******** QMDQT ..... QUOT MIN DEG QUOT TRANSFORM ********/ 8 /***************************************************************/ 9 10 /* PURPOSE - THIS SUBROUTINE PERFORMS THE QUOTIENT GRAPH */ 11 /* TRANSFORMATION AFTER A NODE HAS BEEN ELIMINATED.*/ 12 13 /* INPUT PARAMETERS -*/ 14 /* ROOT - THE NODE JUST ELIMINATED. IT BECOMES THE*/ 15 /* REPRESENTATIVE OF THE NEW SUPERNODE.*/ 16 /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE.*/ 17 /* (RCHSZE, RCHSET) - THE REACHABLE SET OF ROOT IN THE*/ 18 /* OLD QUOTIENT GRAPH.*/ 19 /* NBRHD - THE NEIGHBORHOOD SET WHICH WILL BE MERGED*/ 20 /* WITH ROOT TO FORM THE NEW SUPERNODE.*/ 21 /* MARKER - THE MARKER VECTOR.*/ 22 23 /* UPDATED PARAMETER -*/ 24 /* ADJNCY - BECOMES THE ADJNCY OF THE QUOTIENT GRAPH.*/ 25 /***************************************************************/ 26 PetscErrorCode SPARSEPACKqmdqt(const PetscInt *root, const PetscInt *xadj, const PetscInt *inadjncy, PetscInt *marker, PetscInt *rchsze, PetscInt *rchset, PetscInt *nbrhd) 27 { 28 PetscInt *adjncy = (PetscInt *)inadjncy; /* Used as temporary and reset within this function */ 29 /* System generated locals */ 30 PetscInt i__1, i__2; 31 32 /* Local variables */ 33 PetscInt inhd, irch, node, ilink, j, nabor, jstop, jstrt; 34 35 PetscFunctionBegin; 36 /* Parameter adjustments */ 37 --nbrhd; 38 --rchset; 39 --marker; 40 --adjncy; 41 --xadj; 42 43 irch = 0; 44 inhd = 0; 45 node = *root; 46 L100: 47 jstrt = xadj[node]; 48 jstop = xadj[node + 1] - 2; 49 if (jstop < jstrt) goto L300; 50 51 /* PLACE REACH NODES INTO THE ADJACENT LIST OF NODE*/ 52 i__1 = jstop; 53 for (j = jstrt; j <= i__1; ++j) { 54 ++irch; 55 adjncy[j] = rchset[irch]; 56 if (irch >= *rchsze) goto L400; 57 } 58 /* LINK TO OTHER SPACE PROVIDED BY THE NBRHD SET.*/ 59 L300: 60 ilink = adjncy[jstop + 1]; 61 node = -ilink; 62 if (ilink < 0) goto L100; 63 ++inhd; 64 node = nbrhd[inhd]; 65 adjncy[jstop + 1] = -node; 66 goto L100; 67 /* ALL REACHABLE NODES HAVE BEEN SAVED. END THE ADJ LIST.*/ 68 /* ADD ROOT TO THE NBR LIST OF EACH NODE IN THE REACH SET.*/ 69 L400: 70 adjncy[j + 1] = 0; 71 i__1 = *rchsze; 72 for (irch = 1; irch <= i__1; ++irch) { 73 node = rchset[irch]; 74 if (marker[node] < 0) goto L600; 75 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) goto L500; 82 adjncy[j] = *root; 83 goto L600; 84 L500:; 85 } 86 L600:; 87 } 88 PetscFunctionReturn(PETSC_SUCCESS); 89 } 90