18be712e4SBarry Smith /* genrcm.f -- translated by f2c (version 19931217).*/ 28be712e4SBarry Smith 38be712e4SBarry Smith #include <petscsys.h> 48be712e4SBarry Smith #include <petsc/private/matorderimpl.h> 58be712e4SBarry Smith 68be712e4SBarry Smith /*****************************************************************/ 78be712e4SBarry Smith /*****************************************************************/ 88be712e4SBarry Smith /********* GENRCM ..... GENERAL REVERSE CUTHILL MCKEE ********/ 98be712e4SBarry Smith /*****************************************************************/ 108be712e4SBarry Smith 118be712e4SBarry Smith /* PURPOSE - GENRCM FINDS THE REVERSE CUTHILL-MCKEE*/ 128be712e4SBarry Smith /* ORDERING FOR A GENERAL GRAPH. FOR EACH CONNECTED*/ 138be712e4SBarry Smith /* COMPONENT IN THE GRAPH, GENRCM OBTAINS THE ORDERING*/ 148be712e4SBarry Smith /* BY CALLING THE SUBROUTINE RCM.*/ 158be712e4SBarry Smith 168be712e4SBarry Smith /* INPUT PARAMETERS -*/ 178be712e4SBarry Smith /* NEQNS - NUMBER OF EQUATIONS*/ 188be712e4SBarry Smith /* (XADJ, ADJNCY) - ARRAY PAIR CONTAINING THE ADJACENCY*/ 198be712e4SBarry Smith /* STRUCTURE OF THE GRAPH OF THE MATRIX.*/ 208be712e4SBarry Smith 218be712e4SBarry Smith /* OUTPUT PARAMETER -*/ 228be712e4SBarry Smith /* PERM - VECTOR THAT CONTAINS THE RCM ORDERING.*/ 238be712e4SBarry Smith 248be712e4SBarry Smith /* WORKING PARAMETERS -*/ 258be712e4SBarry Smith /* MASK - IS USED TO MARK VARIABLES THAT HAVE BEEN*/ 268be712e4SBarry Smith /* NUMBERED DURING THE ORDERING PROCESS. IT IS*/ 278be712e4SBarry Smith /* INITIALIZED TO 1, AND SET TO ZERO AS EACH NODE*/ 288be712e4SBarry Smith /* IS NUMBERED.*/ 298be712e4SBarry Smith /* XLS - THE INDEX VECTOR FOR A LEVEL STRUCTURE. THE*/ 308be712e4SBarry Smith /* LEVEL STRUCTURE IS STORED IN THE CURRENTLY*/ 318be712e4SBarry Smith /* UNUSED SPACES IN THE PERMUTATION VECTOR PERM.*/ 328be712e4SBarry Smith 338be712e4SBarry Smith /* PROGRAM SUBROUTINES -*/ 348be712e4SBarry Smith /* FNROOT, RCM.*/ 358be712e4SBarry Smith /*****************************************************************/ 368be712e4SBarry Smith PetscErrorCode SPARSEPACKgenrcm(const PetscInt *neqns, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *perm, PetscInt *mask, PetscInt *xls) 378be712e4SBarry Smith { 388be712e4SBarry Smith /* System generated locals */ 398be712e4SBarry Smith PetscInt i__1; 408be712e4SBarry Smith 418be712e4SBarry Smith /* Local variables */ 428be712e4SBarry Smith PetscInt nlvl, root, i, ccsize; 438be712e4SBarry Smith PetscInt num; 448be712e4SBarry Smith 458be712e4SBarry Smith PetscFunctionBegin; 46*db68c214SStefano Zampini if (!*neqns) PetscFunctionReturn(PETSC_SUCCESS); 47*db68c214SStefano Zampini if (*neqns == 1) { 48*db68c214SStefano Zampini perm[0] = 1; 49*db68c214SStefano Zampini mask[0] = 1; 50*db68c214SStefano Zampini xls[0] = 1; 51*db68c214SStefano Zampini PetscFunctionReturn(PETSC_SUCCESS); 52*db68c214SStefano Zampini } 53*db68c214SStefano Zampini 548be712e4SBarry Smith /* Parameter adjustments */ 558be712e4SBarry Smith --xls; 568be712e4SBarry Smith --mask; 578be712e4SBarry Smith --perm; 588be712e4SBarry Smith --adjncy; 598be712e4SBarry Smith --xadj; 608be712e4SBarry Smith 618be712e4SBarry Smith i__1 = *neqns; 628be712e4SBarry Smith for (i = 1; i <= i__1; ++i) mask[i] = 1; 638be712e4SBarry Smith num = 1; 648be712e4SBarry Smith i__1 = *neqns; 658be712e4SBarry Smith for (i = 1; i <= i__1; ++i) { 668be712e4SBarry Smith /* FOR EACH MASKED CONNECTED COMPONENT ...*/ 678be712e4SBarry Smith if (!mask[i]) goto L200; 688be712e4SBarry Smith root = i; 698be712e4SBarry Smith /* FIRST FIND A PSEUDO-PERIPHERAL NODE ROOT.*/ 708be712e4SBarry Smith /* NOTE THAT THE LEVEL STRUCTURE FOUND BY*/ 718be712e4SBarry Smith /* FNROOT IS STORED STARTING AT PERM(NUM).*/ 728be712e4SBarry Smith /* THEN RCM IS CALLED TO ORDER THE COMPONENT*/ 738be712e4SBarry Smith /* USING ROOT AS THE STARTING NODE.*/ 748be712e4SBarry Smith PetscCall(SPARSEPACKfnroot(&root, &xadj[1], &adjncy[1], &mask[1], &nlvl, &xls[1], &perm[num])); 758be712e4SBarry Smith PetscCall(SPARSEPACKrcm(&root, &xadj[1], &adjncy[1], &mask[1], &perm[num], &ccsize, &xls[1])); 768be712e4SBarry Smith num += ccsize; 778be712e4SBarry Smith if (num > *neqns) PetscFunctionReturn(PETSC_SUCCESS); 788be712e4SBarry Smith L200:; 798be712e4SBarry Smith } 808be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 818be712e4SBarry Smith } 82