1*8be712e4SBarry Smith /* fn1wd.f -- translated by f2c (version 19931217).*/ 2*8be712e4SBarry Smith 3*8be712e4SBarry Smith #include <petsc/private/matorderimpl.h> 4*8be712e4SBarry Smith 5*8be712e4SBarry Smith /*****************************************************************/ 6*8be712e4SBarry Smith /******** FN1WD ..... FIND ONE-WAY DISSECTORS *********/ 7*8be712e4SBarry Smith /*****************************************************************/ 8*8be712e4SBarry Smith /* PURPOSE - THIS SUBROUTINE FINDS ONE-WAY DISSECTORS OF */ 9*8be712e4SBarry Smith /* A CONNECTED COMPONENT SPECIFIED BY MASK AND ROOT. */ 10*8be712e4SBarry Smith /* */ 11*8be712e4SBarry Smith /* INPUT PARAMETERS - */ 12*8be712e4SBarry Smith /* ROOT - A NODE THAT DEFINES (ALONG WITH MASK) THE */ 13*8be712e4SBarry Smith /* COMPONENT TO BE PROCESSED. */ 14*8be712e4SBarry Smith /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE. */ 15*8be712e4SBarry Smith /* */ 16*8be712e4SBarry Smith /* OUTPUT PARAMETERS - */ 17*8be712e4SBarry Smith /* NSEP - NUMBER OF NODES IN THE ONE-WAY DISSECTORS. */ 18*8be712e4SBarry Smith /* SEP - VECTOR CONTAINING THE DISSECTOR NODES. */ 19*8be712e4SBarry Smith /* */ 20*8be712e4SBarry Smith /* UPDATED PARAMETER - */ 21*8be712e4SBarry Smith /* MASK - NODES IN THE DISSECTOR HAVE THEIR MASK VALUES */ 22*8be712e4SBarry Smith /* SET TO ZERO. */ 23*8be712e4SBarry Smith /* */ 24*8be712e4SBarry Smith /* WORKING PARAMETERS- */ 25*8be712e4SBarry Smith /* (XLS, LS) - LEVEL STRUCTURE USED BY THE ROUTINE FNROOT. */ 26*8be712e4SBarry Smith /* */ 27*8be712e4SBarry Smith /* PROGRAM SUBROUTINE - */ 28*8be712e4SBarry Smith /* FNROOT. */ 29*8be712e4SBarry Smith /*****************************************************************/ 30*8be712e4SBarry Smith PetscErrorCode SPARSEPACKfn1wd(PetscInt *root, const PetscInt *inxadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nsep, PetscInt *sep, PetscInt *nlvl, PetscInt *xls, PetscInt *ls) 31*8be712e4SBarry Smith { 32*8be712e4SBarry Smith PetscInt *xadj = (PetscInt *)inxadj; /* Used as temporary and reset */ 33*8be712e4SBarry Smith /* System generated locals */ 34*8be712e4SBarry Smith PetscInt i__1, i__2; 35*8be712e4SBarry Smith 36*8be712e4SBarry Smith /* Local variables */ 37*8be712e4SBarry Smith PetscInt node, i, j, k; 38*8be712e4SBarry Smith PetscReal width, fnlvl; 39*8be712e4SBarry Smith PetscInt kstop, kstrt, lp1beg, lp1end; 40*8be712e4SBarry Smith PetscReal deltp1; 41*8be712e4SBarry Smith PetscInt lvlbeg, lvlend; 42*8be712e4SBarry Smith PetscInt nbr, lvl; 43*8be712e4SBarry Smith 44*8be712e4SBarry Smith PetscFunctionBegin; 45*8be712e4SBarry Smith /* Parameter adjustments */ 46*8be712e4SBarry Smith --ls; 47*8be712e4SBarry Smith --xls; 48*8be712e4SBarry Smith --sep; 49*8be712e4SBarry Smith --mask; 50*8be712e4SBarry Smith --adjncy; 51*8be712e4SBarry Smith --xadj; 52*8be712e4SBarry Smith 53*8be712e4SBarry Smith PetscCall(SPARSEPACKfnroot(root, &xadj[1], &adjncy[1], &mask[1], nlvl, &xls[1], &ls[1])); 54*8be712e4SBarry Smith fnlvl = (PetscReal)(*nlvl); 55*8be712e4SBarry Smith *nsep = xls[*nlvl + 1] - 1; 56*8be712e4SBarry Smith width = (PetscReal)(*nsep) / fnlvl; 57*8be712e4SBarry Smith deltp1 = PetscSqrtReal((width * 3. + 13.) / 2.) + 1.; 58*8be712e4SBarry Smith if (*nsep >= 50 && deltp1 <= fnlvl * .5f) goto L300; 59*8be712e4SBarry Smith 60*8be712e4SBarry Smith /* THE COMPONENT IS TOO SMALL, OR THE LEVEL STRUCTURE */ 61*8be712e4SBarry Smith /* IS VERY LONG AND NARROW. RETURN THE WHOLE COMPONENT.*/ 62*8be712e4SBarry Smith i__1 = *nsep; 63*8be712e4SBarry Smith for (i = 1; i <= i__1; ++i) { 64*8be712e4SBarry Smith node = ls[i]; 65*8be712e4SBarry Smith sep[i] = node; 66*8be712e4SBarry Smith mask[node] = 0; 67*8be712e4SBarry Smith } 68*8be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 69*8be712e4SBarry Smith /* FIND THE PARALLEL DISSECTORS.*/ 70*8be712e4SBarry Smith L300: 71*8be712e4SBarry Smith *nsep = 0; 72*8be712e4SBarry Smith i = 0; 73*8be712e4SBarry Smith L400: 74*8be712e4SBarry Smith ++i; 75*8be712e4SBarry Smith lvl = (PetscInt)((PetscReal)i * deltp1 + .5f); 76*8be712e4SBarry Smith if (lvl >= *nlvl) PetscFunctionReturn(PETSC_SUCCESS); 77*8be712e4SBarry Smith lvlbeg = xls[lvl]; 78*8be712e4SBarry Smith lp1beg = xls[lvl + 1]; 79*8be712e4SBarry Smith lvlend = lp1beg - 1; 80*8be712e4SBarry Smith lp1end = xls[lvl + 2] - 1; 81*8be712e4SBarry Smith i__1 = lp1end; 82*8be712e4SBarry Smith for (j = lp1beg; j <= i__1; ++j) { 83*8be712e4SBarry Smith node = ls[j]; 84*8be712e4SBarry Smith xadj[node] = -xadj[node]; 85*8be712e4SBarry Smith } 86*8be712e4SBarry Smith /* NODES IN LEVEL LVL ARE CHOSEN TO FORM DISSECTOR. */ 87*8be712e4SBarry Smith /* INCLUDE ONLY THOSE WITH NEIGHBORS IN LVL+1 LEVEL. */ 88*8be712e4SBarry Smith /* XADJ IS USED TEMPORARILY TO MARK NODES IN LVL+1. */ 89*8be712e4SBarry Smith i__1 = lvlend; 90*8be712e4SBarry Smith for (j = lvlbeg; j <= i__1; ++j) { 91*8be712e4SBarry Smith node = ls[j]; 92*8be712e4SBarry Smith kstrt = xadj[node]; 93*8be712e4SBarry Smith i__2 = xadj[node + 1]; 94*8be712e4SBarry Smith kstop = (PetscInt)PetscAbsInt(i__2) - 1; 95*8be712e4SBarry Smith i__2 = kstop; 96*8be712e4SBarry Smith for (k = kstrt; k <= i__2; ++k) { 97*8be712e4SBarry Smith nbr = adjncy[k]; 98*8be712e4SBarry Smith if (xadj[nbr] > 0) goto L600; 99*8be712e4SBarry Smith ++(*nsep); 100*8be712e4SBarry Smith sep[*nsep] = node; 101*8be712e4SBarry Smith mask[node] = 0; 102*8be712e4SBarry Smith goto L700; 103*8be712e4SBarry Smith L600:; 104*8be712e4SBarry Smith } 105*8be712e4SBarry Smith L700:; 106*8be712e4SBarry Smith } 107*8be712e4SBarry Smith i__1 = lp1end; 108*8be712e4SBarry Smith for (j = lp1beg; j <= i__1; ++j) { 109*8be712e4SBarry Smith node = ls[j]; 110*8be712e4SBarry Smith xadj[node] = -xadj[node]; 111*8be712e4SBarry Smith } 112*8be712e4SBarry Smith goto L400; 113*8be712e4SBarry Smith } 114