Actual source code: gen1wd.c

  2: /* gen1wd.f -- translated by f2c (version 19931217).*/

  4: #include <petscsys.h>


  9: /*****************************************************************/
 10: /***********     GEN1WD ..... GENERAL ONE-WAY DISSECTION  ********/
 11: /*****************************************************************/

 13: /*    PURPOSE - GEN1WD FINDS A ONE-WAY DISSECTION PARTITIONING*/
 14: /*       FOR A GENERAL GRAPH.  FN1WD IS USED FOR EACH CONNECTED*/
 15: /*       COMPONENT.*/

 17: /*    INPUT PARAMETERS -*/
 18: /*       NEQNS - NUMBER OF EQUATIONS.*/
 19: /*       (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE PAIR.*/

 21: /*    OUTPUT PARAMETERS -*/
 22: /*       (NBLKS, XBLK) - THE PARTITIONING FOUND.*/
 23: /*       PERM - THE ONE-WAY DISSECTION ORDERING.*/

 25: /*    WORKING VECTORS -*/
 26: /*       MASK - IS USED TO MARK VARIABLES THAT HAVE*/
 27: /*              BEEN NUMBERED DURING THE ORDERING PROCESS.*/
 28: /*       (XLS, LS) - LEVEL STRUCTURE USED BY ../../..LS.*/

 30: /*    PROGRAM SUBROUTINES -*/
 31: /*       FN1WD, REVRSE, ../../..LS.*/
 32: /****************************************************************/
 35: PetscErrorCode SPARSEPACKgen1wd(PetscInt *neqns, PetscInt *xadj, PetscInt *adjncy, 
 36:                                 PetscInt *mask, PetscInt *nblks, PetscInt *xblk, PetscInt *perm, PetscInt *xls, PetscInt *ls)
 37: {
 38:     /* System generated locals */
 39:     PetscInt i__1, i__2, i__3;

 41:     /* Local variables */
 42:     PetscInt node, nsep, lnum, nlvl, root;
 43:     PetscInt i, j, k, ccsize;
 44:     PetscInt num;

 47:     /* Parameter adjustments */
 48:     --ls;
 49:     --xls;
 50:     --perm;
 51:     --xblk;
 52:     --mask;
 53:     --xadj;
 54:     --adjncy;

 56:     i__1 = *neqns;
 57:     for (i = 1; i <= i__1; ++i) {
 58:         mask[i] = 1;
 59:     }
 60:     *nblks = 0;
 61:     num = 0;
 62:     i__1 = *neqns;
 63:     for (i = 1; i <= i__1; ++i) {
 64:         if (!mask[i]) {
 65:             goto L400;
 66:         }
 67: /*             FIND A ONE-WAY DISSECTOR FOR EACH COMPONENT.*/
 68:         root = i;
 69:         SPARSEPACKfn1wd(&root, &xadj[1], &adjncy[1], &mask[1], &nsep, &perm[num + 1], &
 70:                 nlvl, &xls[1], &ls[1]);
 71:         num += nsep;
 72:         ++(*nblks);
 73:         xblk[*nblks] = *neqns - num + 1;
 74:         ccsize = xls[nlvl + 1] - 1;
 75: /*             NUMBER THE REMAINING NODES IN THE COMPONENT.*/
 76: /*             EACH COMPONENT IN THE REMAINING SUBGRAPH FORMS*/
 77: /*             A NEW BLOCK IN THE PARTITIONING.*/
 78:         i__2 = ccsize;
 79:         for (j = 1; j <= i__2; ++j) {
 80:             node = ls[j];
 81:             if (!mask[node]) {
 82:                 goto L300;
 83:             }
 84:             SPARSEPACKrootls(&node, &xadj[1], &adjncy[1], &mask[1], &nlvl, &xls[1], &
 85:                     perm[num + 1]);
 86:             lnum = num + 1;
 87:             num = num + xls[nlvl + 1] - 1;
 88:             ++(*nblks);
 89:             xblk[*nblks] = *neqns - num + 1;
 90:             i__3 = num;
 91:             for (k = lnum; k <= i__3; ++k) {
 92:                 node = perm[k];
 93:                 mask[node] = 0;
 94:             }
 95:             if (num > *neqns) {
 96:                 goto L500;
 97:             }
 98: L300:
 99:             ;
100:         }
101: L400:
102:         ;
103:     }
104: /*       SINCE DISSECTORS FOUND FIRST SHOULD BE ORDERED LAST,*/
105: /*       ROUTINE REVRSE IS CALLED TO ADJUST THE ORDERING*/
106: /*       VECTOR, AND THE BLOCK INDEX VECTOR.*/
107: L500:
108:     SPARSEPACKrevrse(neqns, &perm[1]);
109:     SPARSEPACKrevrse(nblks, &xblk[1]);
110:     xblk[*nblks + 1] = *neqns + 1;
111:     return(0);
112: }