1*8be712e4SBarry Smith /* fndsep.f -- translated by f2c (version 19931217). 2*8be712e4SBarry Smith */ 3*8be712e4SBarry Smith 4*8be712e4SBarry Smith #include <petsc/private/matorderimpl.h> 5*8be712e4SBarry Smith 6*8be712e4SBarry Smith /*****************************************************************/ 7*8be712e4SBarry Smith /************* FNDSEP ..... FIND SEPARATOR *************/ 8*8be712e4SBarry Smith /*****************************************************************/ 9*8be712e4SBarry Smith /* PURPOSE - THIS ROUTINE IS USED TO FIND A SMALL */ 10*8be712e4SBarry Smith /* SEPARATOR FOR A CONNECTED COMPONENT SPECIFIED */ 11*8be712e4SBarry Smith /* BY MASK IN THE GIVEN GRAPH. */ 12*8be712e4SBarry Smith /* */ 13*8be712e4SBarry Smith /* INPUT PARAMETERS - */ 14*8be712e4SBarry Smith /* ROOT - IS THE NODE THAT DETERMINES THE MASKED */ 15*8be712e4SBarry Smith /* COMPONENT. */ 16*8be712e4SBarry Smith /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE PAIR. */ 17*8be712e4SBarry Smith /* */ 18*8be712e4SBarry Smith /* OUTPUT PARAMETERS - */ 19*8be712e4SBarry Smith /* NSEP - NUMBER OF VARIABLES IN THE SEPARATOR. */ 20*8be712e4SBarry Smith /* SEP - VECTOR CONTAINING THE SEPARATOR NODES. */ 21*8be712e4SBarry Smith /* */ 22*8be712e4SBarry Smith /* UPDATED PARAMETER - */ 23*8be712e4SBarry Smith /* MASK - NODES IN THE SEPARATOR HAVE THEIR MASK */ 24*8be712e4SBarry Smith /* VALUES SET TO ZERO. */ 25*8be712e4SBarry Smith /* */ 26*8be712e4SBarry Smith /* WORKING PARAMETERS - */ 27*8be712e4SBarry Smith /* (XLS, LS) - LEVEL STRUCTURE PAIR FOR LEVEL STRUCTURE */ 28*8be712e4SBarry Smith /* FOUND BY FNROOT. */ 29*8be712e4SBarry Smith /* */ 30*8be712e4SBarry Smith /* PROGRAM SUBROUTINES - */ 31*8be712e4SBarry Smith /* FNROOT. */ 32*8be712e4SBarry Smith /* */ 33*8be712e4SBarry Smith /*****************************************************************/ 34*8be712e4SBarry Smith PetscErrorCode SPARSEPACKfndsep(PetscInt *root, const PetscInt *inxadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nsep, PetscInt *sep, PetscInt *xls, PetscInt *ls) 35*8be712e4SBarry Smith { 36*8be712e4SBarry Smith PetscInt *xadj = (PetscInt *)inxadj; /* Used as temporary and reset within this function */ 37*8be712e4SBarry Smith 38*8be712e4SBarry Smith /* System generated locals */ 39*8be712e4SBarry Smith PetscInt i__1, i__2; 40*8be712e4SBarry Smith 41*8be712e4SBarry Smith /* Local variables */ 42*8be712e4SBarry Smith PetscInt node, nlvl, i, j, jstop, jstrt, mp1beg, mp1end, midbeg, midend, midlvl; 43*8be712e4SBarry Smith PetscInt nbr; 44*8be712e4SBarry Smith 45*8be712e4SBarry Smith PetscFunctionBegin; 46*8be712e4SBarry Smith /* Parameter adjustments */ 47*8be712e4SBarry Smith --ls; 48*8be712e4SBarry Smith --xls; 49*8be712e4SBarry Smith --sep; 50*8be712e4SBarry Smith --mask; 51*8be712e4SBarry Smith --adjncy; 52*8be712e4SBarry Smith --xadj; 53*8be712e4SBarry Smith 54*8be712e4SBarry Smith PetscCall(SPARSEPACKfnroot(root, &xadj[1], &adjncy[1], &mask[1], &nlvl, &xls[1], &ls[1])); 55*8be712e4SBarry Smith /* IF THE NUMBER OF LEVELS IS LESS THAN 3, RETURN */ 56*8be712e4SBarry Smith /* THE WHOLE COMPONENT AS THE SEPARATOR.*/ 57*8be712e4SBarry Smith if (nlvl >= 3) goto L200; 58*8be712e4SBarry Smith *nsep = xls[nlvl + 1] - 1; 59*8be712e4SBarry Smith i__1 = *nsep; 60*8be712e4SBarry Smith for (i = 1; i <= i__1; ++i) { 61*8be712e4SBarry Smith node = ls[i]; 62*8be712e4SBarry Smith sep[i] = node; 63*8be712e4SBarry Smith mask[node] = 0; 64*8be712e4SBarry Smith } 65*8be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 66*8be712e4SBarry Smith /* FIND THE MIDDLE LEVEL OF THE ROOTED LEVEL STRUCTURE.*/ 67*8be712e4SBarry Smith L200: 68*8be712e4SBarry Smith midlvl = (nlvl + 2) / 2; 69*8be712e4SBarry Smith midbeg = xls[midlvl]; 70*8be712e4SBarry Smith mp1beg = xls[midlvl + 1]; 71*8be712e4SBarry Smith midend = mp1beg - 1; 72*8be712e4SBarry Smith mp1end = xls[midlvl + 2] - 1; 73*8be712e4SBarry Smith /* THE SEPARATOR IS OBTAINED BY INCLUDING ONLY THOSE*/ 74*8be712e4SBarry Smith /* MIDDLE-LEVEL NODES WITH NEIGHBORS IN THE MIDDLE+1*/ 75*8be712e4SBarry Smith /* LEVEL. XADJ IS USED TEMPORARILY TO MARK THOSE*/ 76*8be712e4SBarry Smith /* NODES IN THE MIDDLE+1 LEVEL.*/ 77*8be712e4SBarry Smith i__1 = mp1end; 78*8be712e4SBarry Smith for (i = mp1beg; i <= i__1; ++i) { 79*8be712e4SBarry Smith node = ls[i]; 80*8be712e4SBarry Smith xadj[node] = -xadj[node]; 81*8be712e4SBarry Smith } 82*8be712e4SBarry Smith *nsep = 0; 83*8be712e4SBarry Smith i__1 = midend; 84*8be712e4SBarry Smith for (i = midbeg; i <= i__1; ++i) { 85*8be712e4SBarry Smith node = ls[i]; 86*8be712e4SBarry Smith jstrt = xadj[node]; 87*8be712e4SBarry Smith i__2 = xadj[node + 1]; 88*8be712e4SBarry Smith jstop = (PetscInt)PetscAbsInt(i__2) - 1; 89*8be712e4SBarry Smith i__2 = jstop; 90*8be712e4SBarry Smith for (j = jstrt; j <= i__2; ++j) { 91*8be712e4SBarry Smith nbr = adjncy[j]; 92*8be712e4SBarry Smith if (xadj[nbr] > 0) goto L400; 93*8be712e4SBarry Smith ++(*nsep); 94*8be712e4SBarry Smith sep[*nsep] = node; 95*8be712e4SBarry Smith mask[node] = 0; 96*8be712e4SBarry Smith goto L500; 97*8be712e4SBarry Smith L400:; 98*8be712e4SBarry Smith } 99*8be712e4SBarry Smith L500:; 100*8be712e4SBarry Smith } 101*8be712e4SBarry Smith /* RESET XADJ TO ITS CORRECT SIGN.*/ 102*8be712e4SBarry Smith i__1 = mp1end; 103*8be712e4SBarry Smith for (i = mp1beg; i <= i__1; ++i) { 104*8be712e4SBarry Smith node = ls[i]; 105*8be712e4SBarry Smith xadj[node] = -xadj[node]; 106*8be712e4SBarry Smith } 107*8be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 108*8be712e4SBarry Smith } 109