1*8be712e4SBarry Smith /* gen1wd.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 /*********** GEN1WD ..... GENERAL ONE-WAY DISSECTION ********/ 8*8be712e4SBarry Smith /*****************************************************************/ 9*8be712e4SBarry Smith 10*8be712e4SBarry Smith /* PURPOSE - GEN1WD FINDS A ONE-WAY DISSECTION PARTITIONING*/ 11*8be712e4SBarry Smith /* FOR A GENERAL GRAPH. FN1WD IS USED FOR EACH CONNECTED*/ 12*8be712e4SBarry Smith /* COMPONENT.*/ 13*8be712e4SBarry Smith 14*8be712e4SBarry Smith /* INPUT PARAMETERS -*/ 15*8be712e4SBarry Smith /* NEQNS - NUMBER OF EQUATIONS.*/ 16*8be712e4SBarry Smith /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE PAIR.*/ 17*8be712e4SBarry Smith 18*8be712e4SBarry Smith /* OUTPUT PARAMETERS -*/ 19*8be712e4SBarry Smith /* (NBLKS, XBLK) - THE PARTITIONING FOUND.*/ 20*8be712e4SBarry Smith /* PERM - THE ONE-WAY DISSECTION ORDERING.*/ 21*8be712e4SBarry Smith 22*8be712e4SBarry Smith /* WORKING VECTORS -*/ 23*8be712e4SBarry Smith /* MASK - IS USED TO MARK VARIABLES THAT HAVE*/ 24*8be712e4SBarry Smith /* BEEN NUMBERED DURING THE ORDERING PROCESS.*/ 25*8be712e4SBarry Smith /* (XLS, LS) - LEVEL STRUCTURE USED BY ROOTLS.*/ 26*8be712e4SBarry Smith 27*8be712e4SBarry Smith /* PROGRAM SUBROUTINES -*/ 28*8be712e4SBarry Smith /* FN1WD, REVRSE, ROOTLS.*/ 29*8be712e4SBarry Smith /****************************************************************/ 30*8be712e4SBarry Smith PetscErrorCode SPARSEPACKgen1wd(const PetscInt *neqns, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nblks, PetscInt *xblk, PetscInt *perm, PetscInt *xls, PetscInt *ls) 31*8be712e4SBarry Smith { 32*8be712e4SBarry Smith /* System generated locals */ 33*8be712e4SBarry Smith PetscInt i__1, i__2, i__3; 34*8be712e4SBarry Smith 35*8be712e4SBarry Smith /* Local variables */ 36*8be712e4SBarry Smith PetscInt node, nsep, lnum, nlvl, root; 37*8be712e4SBarry Smith PetscInt i, j, k, ccsize; 38*8be712e4SBarry Smith PetscInt num; 39*8be712e4SBarry Smith 40*8be712e4SBarry Smith PetscFunctionBegin; 41*8be712e4SBarry Smith /* Parameter adjustments */ 42*8be712e4SBarry Smith --ls; 43*8be712e4SBarry Smith --xls; 44*8be712e4SBarry Smith --perm; 45*8be712e4SBarry Smith --xblk; 46*8be712e4SBarry Smith --mask; 47*8be712e4SBarry Smith --xadj; 48*8be712e4SBarry Smith --adjncy; 49*8be712e4SBarry Smith 50*8be712e4SBarry Smith i__1 = *neqns; 51*8be712e4SBarry Smith for (i = 1; i <= i__1; ++i) mask[i] = 1; 52*8be712e4SBarry Smith *nblks = 0; 53*8be712e4SBarry Smith num = 0; 54*8be712e4SBarry Smith i__1 = *neqns; 55*8be712e4SBarry Smith for (i = 1; i <= i__1; ++i) { 56*8be712e4SBarry Smith if (!mask[i]) goto L400; 57*8be712e4SBarry Smith /* FIND A ONE-WAY DISSECTOR FOR EACH COMPONENT.*/ 58*8be712e4SBarry Smith root = i; 59*8be712e4SBarry Smith PetscCall(SPARSEPACKfn1wd(&root, &xadj[1], &adjncy[1], &mask[1], &nsep, &perm[num + 1], &nlvl, &xls[1], &ls[1])); 60*8be712e4SBarry Smith num += nsep; 61*8be712e4SBarry Smith ++(*nblks); 62*8be712e4SBarry Smith xblk[*nblks] = *neqns - num + 1; 63*8be712e4SBarry Smith ccsize = xls[nlvl + 1] - 1; 64*8be712e4SBarry Smith /* NUMBER THE REMAINING NODES IN THE COMPONENT.*/ 65*8be712e4SBarry Smith /* EACH COMPONENT IN THE REMAINING SUBGRAPH FORMS*/ 66*8be712e4SBarry Smith /* A NEW BLOCK IN THE PARTITIONING.*/ 67*8be712e4SBarry Smith i__2 = ccsize; 68*8be712e4SBarry Smith for (j = 1; j <= i__2; ++j) { 69*8be712e4SBarry Smith node = ls[j]; 70*8be712e4SBarry Smith if (!mask[node]) goto L300; 71*8be712e4SBarry Smith PetscCall(SPARSEPACKrootls(&node, &xadj[1], &adjncy[1], &mask[1], &nlvl, &xls[1], &perm[num + 1])); 72*8be712e4SBarry Smith lnum = num + 1; 73*8be712e4SBarry Smith num = num + xls[nlvl + 1] - 1; 74*8be712e4SBarry Smith ++(*nblks); 75*8be712e4SBarry Smith xblk[*nblks] = *neqns - num + 1; 76*8be712e4SBarry Smith i__3 = num; 77*8be712e4SBarry Smith for (k = lnum; k <= i__3; ++k) { 78*8be712e4SBarry Smith node = perm[k]; 79*8be712e4SBarry Smith mask[node] = 0; 80*8be712e4SBarry Smith } 81*8be712e4SBarry Smith if (num > *neqns) goto L500; 82*8be712e4SBarry Smith L300:; 83*8be712e4SBarry Smith } 84*8be712e4SBarry Smith L400:; 85*8be712e4SBarry Smith } 86*8be712e4SBarry Smith /* SINCE DISSECTORS FOUND FIRST SHOULD BE ORDERED LAST,*/ 87*8be712e4SBarry Smith /* ROUTINE REVRSE IS CALLED TO ADJUST THE ORDERING*/ 88*8be712e4SBarry Smith /* VECTOR, AND THE BLOCK INDEX VECTOR.*/ 89*8be712e4SBarry Smith L500: 90*8be712e4SBarry Smith PetscCall(SPARSEPACKrevrse(neqns, &perm[1])); 91*8be712e4SBarry Smith PetscCall(SPARSEPACKrevrse(nblks, &xblk[1])); 92*8be712e4SBarry Smith xblk[*nblks + 1] = *neqns + 1; 93*8be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 94*8be712e4SBarry Smith } 95