Actual source code: gen1wd.c

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

  3: #include <petscsys.h>
  4: #include <petsc/private/matorderimpl.h>

  6: /*****************************************************************/
  7: /***********     GEN1WD ..... GENERAL ONE-WAY DISSECTION  ********/
  8: /*****************************************************************/

 10: /*    PURPOSE - GEN1WD FINDS A ONE-WAY DISSECTION PARTITIONING*/
 11: /*       FOR A GENERAL GRAPH.  FN1WD IS USED FOR EACH CONNECTED*/
 12: /*       COMPONENT.*/

 14: /*    INPUT PARAMETERS -*/
 15: /*       NEQNS - NUMBER OF EQUATIONS.*/
 16: /*       (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE PAIR.*/

 18: /*    OUTPUT PARAMETERS -*/
 19: /*       (NBLKS, XBLK) - THE PARTITIONING FOUND.*/
 20: /*       PERM - THE ONE-WAY DISSECTION ORDERING.*/

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

 27: /*    PROGRAM SUBROUTINES -*/
 28: /*       FN1WD, REVRSE, ROOTLS.*/
 29: /****************************************************************/
 30: PetscErrorCode SPARSEPACKgen1wd(const PetscInt *neqns, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nblks, PetscInt *xblk, PetscInt *perm, PetscInt *xls, PetscInt *ls)
 31: {
 32:   /* System generated locals */
 33:   PetscInt i__1, i__2, i__3;

 35:   /* Local variables */
 36:   PetscInt node, nsep, lnum, nlvl, root;
 37:   PetscInt i, j, k, ccsize;
 38:   PetscInt num;

 40:   PetscFunctionBegin;
 41:   /* Parameter adjustments */
 42:   --ls;
 43:   --xls;
 44:   --perm;
 45:   --xblk;
 46:   --mask;
 47:   --xadj;
 48:   --adjncy;

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