18be712e4SBarry Smith /* fndsep.f -- translated by f2c (version 19931217). 28be712e4SBarry Smith */ 38be712e4SBarry Smith 48be712e4SBarry Smith #include <petsc/private/matorderimpl.h> 58be712e4SBarry Smith 68be712e4SBarry Smith /*****************************************************************/ 78be712e4SBarry Smith /************* FNDSEP ..... FIND SEPARATOR *************/ 88be712e4SBarry Smith /*****************************************************************/ 98be712e4SBarry Smith /* PURPOSE - THIS ROUTINE IS USED TO FIND A SMALL */ 108be712e4SBarry Smith /* SEPARATOR FOR A CONNECTED COMPONENT SPECIFIED */ 118be712e4SBarry Smith /* BY MASK IN THE GIVEN GRAPH. */ 128be712e4SBarry Smith /* */ 138be712e4SBarry Smith /* INPUT PARAMETERS - */ 148be712e4SBarry Smith /* ROOT - IS THE NODE THAT DETERMINES THE MASKED */ 158be712e4SBarry Smith /* COMPONENT. */ 168be712e4SBarry Smith /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE PAIR. */ 178be712e4SBarry Smith /* */ 188be712e4SBarry Smith /* OUTPUT PARAMETERS - */ 198be712e4SBarry Smith /* NSEP - NUMBER OF VARIABLES IN THE SEPARATOR. */ 208be712e4SBarry Smith /* SEP - VECTOR CONTAINING THE SEPARATOR NODES. */ 218be712e4SBarry Smith /* */ 228be712e4SBarry Smith /* UPDATED PARAMETER - */ 238be712e4SBarry Smith /* MASK - NODES IN THE SEPARATOR HAVE THEIR MASK */ 248be712e4SBarry Smith /* VALUES SET TO ZERO. */ 258be712e4SBarry Smith /* */ 268be712e4SBarry Smith /* WORKING PARAMETERS - */ 278be712e4SBarry Smith /* (XLS, LS) - LEVEL STRUCTURE PAIR FOR LEVEL STRUCTURE */ 288be712e4SBarry Smith /* FOUND BY FNROOT. */ 298be712e4SBarry Smith /* */ 308be712e4SBarry Smith /* PROGRAM SUBROUTINES - */ 318be712e4SBarry Smith /* FNROOT. */ 328be712e4SBarry Smith /* */ 338be712e4SBarry Smith /*****************************************************************/ 348be712e4SBarry Smith PetscErrorCode SPARSEPACKfndsep(PetscInt *root, const PetscInt *inxadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nsep, PetscInt *sep, PetscInt *xls, PetscInt *ls) 358be712e4SBarry Smith { 368be712e4SBarry Smith PetscInt *xadj = (PetscInt *)inxadj; /* Used as temporary and reset within this function */ 378be712e4SBarry Smith 388be712e4SBarry Smith /* System generated locals */ 398be712e4SBarry Smith PetscInt i__1, i__2; 408be712e4SBarry Smith 418be712e4SBarry Smith /* Local variables */ 428be712e4SBarry Smith PetscInt node, nlvl, i, j, jstop, jstrt, mp1beg, mp1end, midbeg, midend, midlvl; 438be712e4SBarry Smith PetscInt nbr; 448be712e4SBarry Smith 458be712e4SBarry Smith PetscFunctionBegin; 468be712e4SBarry Smith /* Parameter adjustments */ 478be712e4SBarry Smith --ls; 488be712e4SBarry Smith --xls; 498be712e4SBarry Smith --sep; 508be712e4SBarry Smith --mask; 518be712e4SBarry Smith --adjncy; 528be712e4SBarry Smith --xadj; 538be712e4SBarry Smith 548be712e4SBarry Smith PetscCall(SPARSEPACKfnroot(root, &xadj[1], &adjncy[1], &mask[1], &nlvl, &xls[1], &ls[1])); 558be712e4SBarry Smith /* IF THE NUMBER OF LEVELS IS LESS THAN 3, RETURN */ 568be712e4SBarry Smith /* THE WHOLE COMPONENT AS THE SEPARATOR.*/ 578be712e4SBarry Smith if (nlvl >= 3) goto L200; 588be712e4SBarry Smith *nsep = xls[nlvl + 1] - 1; 598be712e4SBarry Smith i__1 = *nsep; 608be712e4SBarry Smith for (i = 1; i <= i__1; ++i) { 618be712e4SBarry Smith node = ls[i]; 628be712e4SBarry Smith sep[i] = node; 638be712e4SBarry Smith mask[node] = 0; 648be712e4SBarry Smith } 658be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 668be712e4SBarry Smith /* FIND THE MIDDLE LEVEL OF THE ROOTED LEVEL STRUCTURE.*/ 678be712e4SBarry Smith L200: 688be712e4SBarry Smith midlvl = (nlvl + 2) / 2; 698be712e4SBarry Smith midbeg = xls[midlvl]; 708be712e4SBarry Smith mp1beg = xls[midlvl + 1]; 718be712e4SBarry Smith midend = mp1beg - 1; 728be712e4SBarry Smith mp1end = xls[midlvl + 2] - 1; 738be712e4SBarry Smith /* THE SEPARATOR IS OBTAINED BY INCLUDING ONLY THOSE*/ 748be712e4SBarry Smith /* MIDDLE-LEVEL NODES WITH NEIGHBORS IN THE MIDDLE+1*/ 758be712e4SBarry Smith /* LEVEL. XADJ IS USED TEMPORARILY TO MARK THOSE*/ 768be712e4SBarry Smith /* NODES IN THE MIDDLE+1 LEVEL.*/ 778be712e4SBarry Smith i__1 = mp1end; 788be712e4SBarry Smith for (i = mp1beg; i <= i__1; ++i) { 798be712e4SBarry Smith node = ls[i]; 808be712e4SBarry Smith xadj[node] = -xadj[node]; 818be712e4SBarry Smith } 828be712e4SBarry Smith *nsep = 0; 838be712e4SBarry Smith i__1 = midend; 848be712e4SBarry Smith for (i = midbeg; i <= i__1; ++i) { 858be712e4SBarry Smith node = ls[i]; 868be712e4SBarry Smith jstrt = xadj[node]; 878be712e4SBarry Smith i__2 = xadj[node + 1]; 88*835f2295SStefano Zampini jstop = PetscAbsInt(i__2) - 1; 898be712e4SBarry Smith i__2 = jstop; 908be712e4SBarry Smith for (j = jstrt; j <= i__2; ++j) { 918be712e4SBarry Smith nbr = adjncy[j]; 928be712e4SBarry Smith if (xadj[nbr] > 0) goto L400; 938be712e4SBarry Smith ++(*nsep); 948be712e4SBarry Smith sep[*nsep] = node; 958be712e4SBarry Smith mask[node] = 0; 968be712e4SBarry Smith goto L500; 978be712e4SBarry Smith L400:; 988be712e4SBarry Smith } 998be712e4SBarry Smith L500:; 1008be712e4SBarry Smith } 1018be712e4SBarry Smith /* RESET XADJ TO ITS CORRECT SIGN.*/ 1028be712e4SBarry Smith i__1 = mp1end; 1038be712e4SBarry Smith for (i = mp1beg; i <= i__1; ++i) { 1048be712e4SBarry Smith node = ls[i]; 1058be712e4SBarry Smith xadj[node] = -xadj[node]; 1068be712e4SBarry Smith } 1078be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 1088be712e4SBarry Smith } 109