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