Actual source code: genrcm.c

  1: /* genrcm.f -- translated by f2c (version 19931217).*/

  3: #include <petscsys.h>
  4: #include <petsc/private/matorderimpl.h>

  6: /*****************************************************************/
  7: /*****************************************************************/
  8: /*********   GENRCM ..... GENERAL REVERSE CUTHILL MCKEE   ********/
  9: /*****************************************************************/

 11: /*    PURPOSE - GENRCM FINDS THE REVERSE CUTHILL-MCKEE*/
 12: /*       ORDERING FOR A GENERAL GRAPH. FOR EACH CONNECTED*/
 13: /*       COMPONENT IN THE GRAPH, GENRCM OBTAINS THE ORDERING*/
 14: /*       BY CALLING THE SUBROUTINE RCM.*/

 16: /*    INPUT PARAMETERS -*/
 17: /*       NEQNS - NUMBER OF EQUATIONS*/
 18: /*       (XADJ, ADJNCY) - ARRAY PAIR CONTAINING THE ADJACENCY*/
 19: /*              STRUCTURE OF THE GRAPH OF THE MATRIX.*/

 21: /*    OUTPUT PARAMETER -*/
 22: /*       PERM - VECTOR THAT CONTAINS THE RCM ORDERING.*/

 24: /*    WORKING PARAMETERS -*/
 25: /*       MASK - IS USED TO MARK VARIABLES THAT HAVE BEEN*/
 26: /*              NUMBERED DURING THE ORDERING PROCESS. IT IS*/
 27: /*              INITIALIZED TO 1, AND SET TO ZERO AS EACH NODE*/
 28: /*              IS NUMBERED.*/
 29: /*       XLS - THE INDEX VECTOR FOR A LEVEL STRUCTURE.  THE*/
 30: /*              LEVEL STRUCTURE IS STORED IN THE CURRENTLY*/
 31: /*              UNUSED SPACES IN THE PERMUTATION VECTOR PERM.*/

 33: /*    PROGRAM SUBROUTINES -*/
 34: /*       FNROOT, RCM.*/
 35: /*****************************************************************/
 36: PetscErrorCode SPARSEPACKgenrcm(const PetscInt *neqns, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *perm, PetscInt *mask, PetscInt *xls)
 37: {
 38:   /* System generated locals */
 39:   PetscInt i__1;

 41:   /* Local variables */
 42:   PetscInt nlvl, root, i, ccsize;
 43:   PetscInt num;

 45:   PetscFunctionBegin;
 46:   /* Parameter adjustments */
 47:   --xls;
 48:   --mask;
 49:   --perm;
 50:   --adjncy;
 51:   --xadj;

 53:   i__1 = *neqns;
 54:   for (i = 1; i <= i__1; ++i) mask[i] = 1;
 55:   num  = 1;
 56:   i__1 = *neqns;
 57:   for (i = 1; i <= i__1; ++i) {
 58:     /*          FOR EACH MASKED CONNECTED COMPONENT ...*/
 59:     if (!mask[i]) goto L200;
 60:     root = i;
 61:     /*             FIRST FIND A PSEUDO-PERIPHERAL NODE ROOT.*/
 62:     /*             NOTE THAT THE LEVEL STRUCTURE FOUND BY*/
 63:     /*             FNROOT IS STORED STARTING AT PERM(NUM).*/
 64:     /*             THEN RCM IS CALLED TO ORDER THE COMPONENT*/
 65:     /*             USING ROOT AS THE STARTING NODE.*/
 66:     PetscCall(SPARSEPACKfnroot(&root, &xadj[1], &adjncy[1], &mask[1], &nlvl, &xls[1], &perm[num]));
 67:     PetscCall(SPARSEPACKrcm(&root, &xadj[1], &adjncy[1], &mask[1], &perm[num], &ccsize, &xls[1]));
 68:     num += ccsize;
 69:     if (num > *neqns) PetscFunctionReturn(PETSC_SUCCESS);
 70:   L200:;
 71:   }
 72:   PetscFunctionReturn(PETSC_SUCCESS);
 73: }