1*8be712e4SBarry Smith /* rcm.f -- translated by f2c (version 19931217).*/ 2*8be712e4SBarry Smith 3*8be712e4SBarry Smith #include <petscsys.h> 4*8be712e4SBarry Smith #include <petsc/private/matorderimpl.h> 5*8be712e4SBarry Smith 6*8be712e4SBarry Smith /*****************************************************************/ 7*8be712e4SBarry Smith /********* RCM ..... REVERSE CUTHILL-MCKEE ORDERING *******/ 8*8be712e4SBarry Smith /*****************************************************************/ 9*8be712e4SBarry Smith /* PURPOSE - RCM NUMBERS A CONNECTED COMPONENT SPECIFIED BY */ 10*8be712e4SBarry Smith /* MASK AND ROOT, USING THE RCM ALGORITHM. */ 11*8be712e4SBarry Smith /* THE NUMBERING IS TO BE STARTED AT THE NODE ROOT. */ 12*8be712e4SBarry Smith /* */ 13*8be712e4SBarry Smith /* INPUT PARAMETERS - */ 14*8be712e4SBarry Smith /* ROOT - IS THE NODE THAT DEFINES THE CONNECTED */ 15*8be712e4SBarry Smith /* COMPONENT AND IT IS USED AS THE STARTING */ 16*8be712e4SBarry Smith /* NODE FOR THE RCM ORDERING. */ 17*8be712e4SBarry Smith /* (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR FOR */ 18*8be712e4SBarry Smith /* THE GRAPH. */ 19*8be712e4SBarry Smith /* */ 20*8be712e4SBarry Smith /* UPDATED PARAMETERS - */ 21*8be712e4SBarry Smith /* MASK - ONLY THOSE NODES WITH NONZERO INPUT MASK */ 22*8be712e4SBarry Smith /* VALUES ARE CONSIDERED BY THE ROUTINE. THE */ 23*8be712e4SBarry Smith /* NODES NUMBERED BY RCM WILL HAVE THEIR */ 24*8be712e4SBarry Smith /* MASK VALUES SET TO ZERO. */ 25*8be712e4SBarry Smith /* */ 26*8be712e4SBarry Smith /* OUTPUT PARAMETERS - */ 27*8be712e4SBarry Smith /* PERM - WILL CONTAIN THE RCM ORDERING. */ 28*8be712e4SBarry Smith /* CCSIZE - IS THE SIZE OF THE CONNECTED COMPONENT */ 29*8be712e4SBarry Smith /* THAT HAS BEEN NUMBERED BY RCM. */ 30*8be712e4SBarry Smith /* */ 31*8be712e4SBarry Smith /* WORKING PARAMETER - */ 32*8be712e4SBarry Smith /* DEG - IS A TEMPORARY VECTOR USED TO HOLD THE DEGREE */ 33*8be712e4SBarry Smith /* OF THE NODES IN THE SECTION GRAPH SPECIFIED */ 34*8be712e4SBarry Smith /* BY MASK AND ROOT. */ 35*8be712e4SBarry Smith /* */ 36*8be712e4SBarry Smith /* PROGRAM SUBROUTINES - */ 37*8be712e4SBarry Smith /* DEGREE. */ 38*8be712e4SBarry Smith /* */ 39*8be712e4SBarry Smith /****************************************************************/ 40*8be712e4SBarry Smith PetscErrorCode SPARSEPACKrcm(const PetscInt *root, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *perm, PetscInt *ccsize, PetscInt *deg) 41*8be712e4SBarry Smith { 42*8be712e4SBarry Smith /* System generated locals */ 43*8be712e4SBarry Smith PetscInt i__1, i__2; 44*8be712e4SBarry Smith 45*8be712e4SBarry Smith /* Local variables */ 46*8be712e4SBarry Smith PetscInt node, fnbr, lnbr, i, j, k, l, lperm, jstop, jstrt; 47*8be712e4SBarry Smith PetscInt lbegin, lvlend, nbr; 48*8be712e4SBarry Smith 49*8be712e4SBarry Smith /* FIND THE DEGREES OF THE NODES IN THE */ 50*8be712e4SBarry Smith /* COMPONENT SPECIFIED BY MASK AND ROOT. */ 51*8be712e4SBarry Smith /* ------------------------------------- */ 52*8be712e4SBarry Smith 53*8be712e4SBarry Smith PetscFunctionBegin; 54*8be712e4SBarry Smith /* Parameter adjustments */ 55*8be712e4SBarry Smith --deg; 56*8be712e4SBarry Smith --perm; 57*8be712e4SBarry Smith --mask; 58*8be712e4SBarry Smith --adjncy; 59*8be712e4SBarry Smith --xadj; 60*8be712e4SBarry Smith 61*8be712e4SBarry Smith PetscCall(SPARSEPACKdegree(root, &xadj[1], &adjncy[1], &mask[1], °[1], ccsize, &perm[1])); 62*8be712e4SBarry Smith mask[*root] = 0; 63*8be712e4SBarry Smith if (*ccsize <= 1) PetscFunctionReturn(PETSC_SUCCESS); 64*8be712e4SBarry Smith lvlend = 0; 65*8be712e4SBarry Smith lnbr = 1; 66*8be712e4SBarry Smith /* LBEGIN AND LVLEND POINT TO THE BEGINNING AND */ 67*8be712e4SBarry Smith /* THE END OF THE CURRENT LEVEL RESPECTIVELY. */ 68*8be712e4SBarry Smith L100: 69*8be712e4SBarry Smith lbegin = lvlend + 1; 70*8be712e4SBarry Smith lvlend = lnbr; 71*8be712e4SBarry Smith i__1 = lvlend; 72*8be712e4SBarry Smith for (i = lbegin; i <= i__1; ++i) { 73*8be712e4SBarry Smith /* FOR EACH NODE IN CURRENT LEVEL ... */ 74*8be712e4SBarry Smith node = perm[i]; 75*8be712e4SBarry Smith jstrt = xadj[node]; 76*8be712e4SBarry Smith jstop = xadj[node + 1] - 1; 77*8be712e4SBarry Smith 78*8be712e4SBarry Smith /* FIND THE UNNUMBERED NEIGHBORS OF NODE. */ 79*8be712e4SBarry Smith /* FNBR AND LNBR POINT TO THE FIRST AND LAST */ 80*8be712e4SBarry Smith /* UNNUMBERED NEIGHBORS RESPECTIVELY OF THE CURRENT */ 81*8be712e4SBarry Smith /* NODE IN PERM. */ 82*8be712e4SBarry Smith fnbr = lnbr + 1; 83*8be712e4SBarry Smith i__2 = jstop; 84*8be712e4SBarry Smith for (j = jstrt; j <= i__2; ++j) { 85*8be712e4SBarry Smith nbr = adjncy[j]; 86*8be712e4SBarry Smith if (!mask[nbr]) goto L200; 87*8be712e4SBarry Smith ++lnbr; 88*8be712e4SBarry Smith mask[nbr] = 0; 89*8be712e4SBarry Smith perm[lnbr] = nbr; 90*8be712e4SBarry Smith L200:; 91*8be712e4SBarry Smith } 92*8be712e4SBarry Smith if (fnbr >= lnbr) goto L600; 93*8be712e4SBarry Smith 94*8be712e4SBarry Smith /* SORT THE NEIGHBORS OF NODE IN INCREASING */ 95*8be712e4SBarry Smith /* ORDER BY DEGREE. LINEAR INSERTION IS USED.*/ 96*8be712e4SBarry Smith k = fnbr; 97*8be712e4SBarry Smith L300: 98*8be712e4SBarry Smith l = k; 99*8be712e4SBarry Smith ++k; 100*8be712e4SBarry Smith nbr = perm[k]; 101*8be712e4SBarry Smith L400: 102*8be712e4SBarry Smith if (l < fnbr) goto L500; 103*8be712e4SBarry Smith lperm = perm[l]; 104*8be712e4SBarry Smith if (deg[lperm] <= deg[nbr]) goto L500; 105*8be712e4SBarry Smith perm[l + 1] = lperm; 106*8be712e4SBarry Smith --l; 107*8be712e4SBarry Smith goto L400; 108*8be712e4SBarry Smith L500: 109*8be712e4SBarry Smith perm[l + 1] = nbr; 110*8be712e4SBarry Smith if (k < lnbr) goto L300; 111*8be712e4SBarry Smith L600:; 112*8be712e4SBarry Smith } 113*8be712e4SBarry Smith if (lnbr > lvlend) goto L100; 114*8be712e4SBarry Smith 115*8be712e4SBarry Smith /* WE NOW HAVE THE CUTHILL MCKEE ORDERING.*/ 116*8be712e4SBarry Smith /* REVERSE IT BELOW ...*/ 117*8be712e4SBarry Smith k = *ccsize / 2; 118*8be712e4SBarry Smith l = *ccsize; 119*8be712e4SBarry Smith i__1 = k; 120*8be712e4SBarry Smith for (i = 1; i <= i__1; ++i) { 121*8be712e4SBarry Smith lperm = perm[l]; 122*8be712e4SBarry Smith perm[l] = perm[i]; 123*8be712e4SBarry Smith perm[i] = lperm; 124*8be712e4SBarry Smith --l; 125*8be712e4SBarry Smith } 126*8be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 127*8be712e4SBarry Smith } 128