Actual source code: fn1wd.c


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

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

  6: /*****************************************************************/
  7: /********     FN1WD ..... FIND ONE-WAY DISSECTORS        *********/
  8: /*****************************************************************/
  9: /*    PURPOSE - THIS SUBROUTINE FINDS ONE-WAY DISSECTORS OF      */
 10: /*       A CONNECTED COMPONENT SPECIFIED BY MASK AND ../../...       */
 11: /*                                                               */
 12: /*    INPUT PARAMETERS -                                         */
 13: /*       ../../.. - A NODE THAT DEFINES (ALONG WITH MASK) THE        */
 14: /*              COMPONENT TO BE PROCESSED.                       */
 15: /*       (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE.               */
 16: /*                                                               */
 17: /*    OUTPUT PARAMETERS -                                        */
 18: /*       NSEP - NUMBER OF NODES IN THE ONE-WAY DISSECTORS.       */
 19: /*       SEP - VECTOR CONTAINING THE DISSECTOR NODES.            */
 20: /*                                                               */
 21: /*    UPDATED PARAMETER -                                        */
 22: /*       MASK - NODES IN THE DISSECTOR HAVE THEIR MASK VALUES    */
 23: /*              SET TO ZERO.                                     */
 24: /*                                                               */
 25: /*    WORKING PARAMETERS-                                        */
 26: /*       (XLS, LS) - LEVEL STRUCTURE USED BY THE ROUTINE FN../../... */
 27: /*                                                               */
 28: /*    PROGRAM SUBROUTINE -                                       */
 29: /*       FN../../...                                                 */
 30: /*****************************************************************/
 31: PetscErrorCode SPARSEPACKfn1wd(PetscInt *root, const PetscInt *inxadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nsep, PetscInt *sep, PetscInt *nlvl, PetscInt *xls, PetscInt *ls)
 32: {
 33:   PetscInt *xadj = (PetscInt *)inxadj; /* Used as temporary and reset */
 34:   /* System generated locals */
 35:   PetscInt i__1, i__2;

 37:   /* Local variables */
 38:   PetscInt  node, i, j, k;
 39:   PetscReal width, fnlvl;
 40:   PetscInt  kstop, kstrt, lp1beg, lp1end;
 41:   PetscReal deltp1;
 42:   PetscInt  lvlbeg, lvlend;
 43:   PetscInt  nbr, lvl;

 45:   PetscFunctionBegin;
 46:   /* Parameter adjustments */
 47:   --ls;
 48:   --xls;
 49:   --sep;
 50:   --mask;
 51:   --adjncy;
 52:   --xadj;

 54:   PetscCall(SPARSEPACKfnroot(root, &xadj[1], &adjncy[1], &mask[1], nlvl, &xls[1], &ls[1]));
 55:   fnlvl  = (PetscReal)(*nlvl);
 56:   *nsep  = xls[*nlvl + 1] - 1;
 57:   width  = (PetscReal)(*nsep) / fnlvl;
 58:   deltp1 = PetscSqrtReal((width * 3. + 13.) / 2.) + 1.;
 59:   if (*nsep >= 50 && deltp1 <= fnlvl * .5f) goto L300;

 61:   /*       THE COMPONENT IS TOO SMALL, OR THE LEVEL STRUCTURE */
 62:   /*       IS VERY LONG AND NARROW. RETURN THE WHOLE COMPONENT.*/
 63:   i__1 = *nsep;
 64:   for (i = 1; i <= i__1; ++i) {
 65:     node       = ls[i];
 66:     sep[i]     = node;
 67:     mask[node] = 0;
 68:   }
 69:   PetscFunctionReturn(PETSC_SUCCESS);
 70: /*       FIND THE PARALLEL DISSECTORS.*/
 71: L300:
 72:   *nsep = 0;
 73:   i     = 0;
 74: L400:
 75:   ++i;
 76:   lvl = (PetscInt)((PetscReal)i * deltp1 + .5f);
 77:   if (lvl >= *nlvl) PetscFunctionReturn(PETSC_SUCCESS);
 78:   lvlbeg = xls[lvl];
 79:   lp1beg = xls[lvl + 1];
 80:   lvlend = lp1beg - 1;
 81:   lp1end = xls[lvl + 2] - 1;
 82:   i__1   = lp1end;
 83:   for (j = lp1beg; j <= i__1; ++j) {
 84:     node       = ls[j];
 85:     xadj[node] = -xadj[node];
 86:   }
 87:   /*          NODES IN LEVEL LVL ARE CHOSEN TO FORM DISSECTOR. */
 88:   /*          INCLUDE ONLY THOSE WITH NEIGHBORS IN LVL+1 LEVEL. */
 89:   /*          XADJ IS USED TEMPORARILY TO MARK NODES IN LVL+1.  */
 90:   i__1 = lvlend;
 91:   for (j = lvlbeg; j <= i__1; ++j) {
 92:     node  = ls[j];
 93:     kstrt = xadj[node];
 94:     i__2  = xadj[node + 1];
 95:     kstop = (PetscInt)PetscAbsInt(i__2) - 1;
 96:     i__2  = kstop;
 97:     for (k = kstrt; k <= i__2; ++k) {
 98:       nbr = adjncy[k];
 99:       if (xadj[nbr] > 0) goto L600;
100:       ++(*nsep);
101:       sep[*nsep] = node;
102:       mask[node] = 0;
103:       goto L700;
104:     L600:;
105:     }
106:   L700:;
107:   }
108:   i__1 = lp1end;
109:   for (j = lp1beg; j <= i__1; ++j) {
110:     node       = ls[j];
111:     xadj[node] = -xadj[node];
112:   }
113:   goto L400;
114: }