1*8be712e4SBarry Smith /* fnroot.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 /******** FNROOT ..... FIND PSEUDO-PERIPHERAL NODE ********/ 8*8be712e4SBarry Smith /*****************************************************************/ 9*8be712e4SBarry Smith /* PURPOSE - FNROOT IMPLEMENTS A MODIFIED VERSION OF THE */ 10*8be712e4SBarry Smith /* SCHEME BY GIBBS, POOLE, AND STOCKMEYER TO FIND PSEUDO- */ 11*8be712e4SBarry Smith /* PERIPHERAL NODES. IT DETERMINES SUCH A NODE FOR THE */ 12*8be712e4SBarry Smith /* SECTION SUBGRAPH SPECIFIED BY MASK AND ROOT. */ 13*8be712e4SBarry Smith /* INPUT PARAMETERS - */ 14*8be712e4SBarry Smith /* (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR FOR THE GRAPH. */ 15*8be712e4SBarry Smith /* MASK - SPECIFIES A SECTION SUBGRAPH. NODES FOR WHICH */ 16*8be712e4SBarry Smith /* MASK IS ZERO ARE IGNORED BY FNROOT. */ 17*8be712e4SBarry Smith /* UPDATED PARAMETER - */ 18*8be712e4SBarry Smith /* ROOT - ON INPUT, IT (ALONG WITH MASK) DEFINES THE */ 19*8be712e4SBarry Smith /* COMPONENT FOR WHICH A PSEUDO-PERIPHERAL NODE IS */ 20*8be712e4SBarry Smith /* TO BE FOUND. ON OUTPUT, IT IS THE NODE OBTAINED. */ 21*8be712e4SBarry Smith /* */ 22*8be712e4SBarry Smith /* OUTPUT PARAMETERS - */ 23*8be712e4SBarry Smith /* NLVL - IS THE NUMBER OF LEVELS IN THE LEVEL STRUCTURE */ 24*8be712e4SBarry Smith /* ROOTED AT THE NODE ROOT. */ 25*8be712e4SBarry Smith /* (XLS,LS) - THE LEVEL STRUCTURE ARRAY PAIR CONTAINING */ 26*8be712e4SBarry Smith /* THE LEVEL STRUCTURE FOUND. */ 27*8be712e4SBarry Smith /* */ 28*8be712e4SBarry Smith /* PROGRAM SUBROUTINES - */ 29*8be712e4SBarry Smith /* ROOTLS. */ 30*8be712e4SBarry Smith /* */ 31*8be712e4SBarry Smith /****************************************************************/ 32*8be712e4SBarry Smith PetscErrorCode SPARSEPACKfnroot(PetscInt *root, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nlvl, PetscInt *xls, PetscInt *ls) 33*8be712e4SBarry Smith { 34*8be712e4SBarry Smith /* System generated locals */ 35*8be712e4SBarry Smith PetscInt i__1, i__2; 36*8be712e4SBarry Smith 37*8be712e4SBarry Smith /* Local variables */ 38*8be712e4SBarry Smith PetscInt ndeg, node, j, k, nabor, kstop, jstrt, kstrt, mindeg, ccsize, nunlvl; 39*8be712e4SBarry Smith /* DETERMINE THE LEVEL STRUCTURE ROOTED AT ROOT. */ 40*8be712e4SBarry Smith 41*8be712e4SBarry Smith PetscFunctionBegin; 42*8be712e4SBarry Smith /* Parameter adjustments */ 43*8be712e4SBarry Smith --ls; 44*8be712e4SBarry Smith --xls; 45*8be712e4SBarry Smith --mask; 46*8be712e4SBarry Smith --adjncy; 47*8be712e4SBarry Smith --xadj; 48*8be712e4SBarry Smith 49*8be712e4SBarry Smith PetscCall(SPARSEPACKrootls(root, &xadj[1], &adjncy[1], &mask[1], nlvl, &xls[1], &ls[1])); 50*8be712e4SBarry Smith ccsize = xls[*nlvl + 1] - 1; 51*8be712e4SBarry Smith if (*nlvl == 1 || *nlvl == ccsize) PetscFunctionReturn(PETSC_SUCCESS); 52*8be712e4SBarry Smith 53*8be712e4SBarry Smith /* PICK A NODE WITH MINIMUM DEGREE FROM THE LAST LEVEL.*/ 54*8be712e4SBarry Smith L100: 55*8be712e4SBarry Smith jstrt = xls[*nlvl]; 56*8be712e4SBarry Smith mindeg = ccsize; 57*8be712e4SBarry Smith *root = ls[jstrt]; 58*8be712e4SBarry Smith if (ccsize == jstrt) goto L400; 59*8be712e4SBarry Smith i__1 = ccsize; 60*8be712e4SBarry Smith for (j = jstrt; j <= i__1; ++j) { 61*8be712e4SBarry Smith node = ls[j]; 62*8be712e4SBarry Smith ndeg = 0; 63*8be712e4SBarry Smith kstrt = xadj[node]; 64*8be712e4SBarry Smith kstop = xadj[node + 1] - 1; 65*8be712e4SBarry Smith i__2 = kstop; 66*8be712e4SBarry Smith for (k = kstrt; k <= i__2; ++k) { 67*8be712e4SBarry Smith nabor = adjncy[k]; 68*8be712e4SBarry Smith if (mask[nabor] > 0) ++ndeg; 69*8be712e4SBarry Smith } 70*8be712e4SBarry Smith if (ndeg >= mindeg) goto L300; 71*8be712e4SBarry Smith *root = node; 72*8be712e4SBarry Smith mindeg = ndeg; 73*8be712e4SBarry Smith L300:; 74*8be712e4SBarry Smith } 75*8be712e4SBarry Smith /* AND GENERATE ITS ROOTED LEVEL STRUCTURE.*/ 76*8be712e4SBarry Smith L400: 77*8be712e4SBarry Smith PetscCall(SPARSEPACKrootls(root, &xadj[1], &adjncy[1], &mask[1], &nunlvl, &xls[1], &ls[1])); 78*8be712e4SBarry Smith if (nunlvl <= *nlvl) PetscFunctionReturn(PETSC_SUCCESS); 79*8be712e4SBarry Smith *nlvl = nunlvl; 80*8be712e4SBarry Smith if (*nlvl < ccsize) goto L100; 81*8be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 82*8be712e4SBarry Smith } 83