Actual source code: fn1wd.c

  1: #define PETSCMAT_DLL

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

 5:  #include petsc.h
 6:  #include src/mat/order/order.h
  7: EXTERN PetscErrorCode SPARSEPACKfnroot(PetscInt*, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *);

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

 43:     /* Local variables */
 44:     PetscInt  node, i, j, k;
 45:     PetscReal width, fnlvl;
 46:     PetscInt  kstop, kstrt, lp1beg, lp1end;
 47:     PetscReal deltp1;
 48:     PetscInt  lvlbeg, lvlend;
 49:     PetscInt  nbr, lvl;

 52:     /* Parameter adjustments */
 53:     --ls;
 54:     --xls;
 55:     --sep;
 56:     --mask;
 57:     --adjncy;
 58:     --xadj;

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