Actual source code: gen1wd.c

  1: #define PETSCMAT_DLL

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

 5:  #include petsc.h

  7: EXTERN PetscErrorCode SPARSEPACKfn1wd(PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
  8: EXTERN PetscErrorCode SPARSEPACKrevrse(PetscInt*,PetscInt*),SPARSEPACKrootls(PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);

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

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

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

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

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

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

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

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

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