Actual source code: fndsep.c

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

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

  6: /*****************************************************************/
  7: /*************     FNDSEP ..... FIND SEPARATOR       *************/
  8: /*****************************************************************/
  9: /*    PURPOSE - THIS ROUTINE IS USED TO FIND A SMALL             */
 10: /*              SEPARATOR FOR A CONNECTED COMPONENT SPECIFIED    */
 11: /*              BY MASK IN THE GIVEN GRAPH.                      */
 12: /*                                                               */
 13: /*    INPUT PARAMETERS -                                         */
 14: /*       ROOT - IS THE NODE THAT DETERMINES THE MASKED           */
 15: /*              COMPONENT.                                       */
 16: /*       (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE PAIR.          */
 17: /*                                                               */
 18: /*    OUTPUT PARAMETERS -                                        */
 19: /*       NSEP - NUMBER OF VARIABLES IN THE SEPARATOR.            */
 20: /*       SEP - VECTOR CONTAINING THE SEPARATOR NODES.            */
 21: /*                                                               */
 22: /*    UPDATED PARAMETER -                                        */
 23: /*       MASK - NODES IN THE SEPARATOR HAVE THEIR MASK           */
 24: /*              VALUES SET TO ZERO.                              */
 25: /*                                                               */
 26: /*    WORKING PARAMETERS -                                       */
 27: /*       (XLS, LS) - LEVEL STRUCTURE PAIR FOR LEVEL STRUCTURE    */
 28: /*              FOUND BY FNROOT.                                 */
 29: /*                                                               */
 30: /*    PROGRAM SUBROUTINES -                                      */
 31: /*       FNROOT.                                                 */
 32: /*                                                               */
 33: /*****************************************************************/
 34: PetscErrorCode SPARSEPACKfndsep(PetscInt *root, const PetscInt *inxadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nsep, PetscInt *sep, PetscInt *xls, PetscInt *ls)
 35: {
 36:   PetscInt *xadj = (PetscInt *)inxadj; /* Used as temporary and reset within this function */

 38:   /* System generated locals */
 39:   PetscInt i__1, i__2;

 41:   /* Local variables */
 42:   PetscInt node, nlvl, i, j, jstop, jstrt, mp1beg, mp1end, midbeg, midend, midlvl;
 43:   PetscInt nbr;

 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:   /*       IF THE NUMBER OF LEVELS IS LESS THAN 3, RETURN */
 56:   /*       THE WHOLE COMPONENT AS THE SEPARATOR.*/
 57:   if (nlvl >= 3) goto L200;
 58:   *nsep = xls[nlvl + 1] - 1;
 59:   i__1  = *nsep;
 60:   for (i = 1; i <= i__1; ++i) {
 61:     node       = ls[i];
 62:     sep[i]     = node;
 63:     mask[node] = 0;
 64:   }
 65:   PetscFunctionReturn(PETSC_SUCCESS);
 66: /*       FIND THE MIDDLE LEVEL OF THE ROOTED LEVEL STRUCTURE.*/
 67: L200:
 68:   midlvl = (nlvl + 2) / 2;
 69:   midbeg = xls[midlvl];
 70:   mp1beg = xls[midlvl + 1];
 71:   midend = mp1beg - 1;
 72:   mp1end = xls[midlvl + 2] - 1;
 73:   /*       THE SEPARATOR IS OBTAINED BY INCLUDING ONLY THOSE*/
 74:   /*       MIDDLE-LEVEL NODES WITH NEIGHBORS IN THE MIDDLE+1*/
 75:   /*       LEVEL. XADJ IS USED TEMPORARILY TO MARK THOSE*/
 76:   /*       NODES IN THE MIDDLE+1 LEVEL.*/
 77:   i__1 = mp1end;
 78:   for (i = mp1beg; i <= i__1; ++i) {
 79:     node       = ls[i];
 80:     xadj[node] = -xadj[node];
 81:   }
 82:   *nsep = 0;
 83:   i__1  = midend;
 84:   for (i = midbeg; i <= i__1; ++i) {
 85:     node  = ls[i];
 86:     jstrt = xadj[node];
 87:     i__2  = xadj[node + 1];
 88:     jstop = (PetscInt)PetscAbsInt(i__2) - 1;
 89:     i__2  = jstop;
 90:     for (j = jstrt; j <= i__2; ++j) {
 91:       nbr = adjncy[j];
 92:       if (xadj[nbr] > 0) goto L400;
 93:       ++(*nsep);
 94:       sep[*nsep] = node;
 95:       mask[node] = 0;
 96:       goto L500;
 97:     L400:;
 98:     }
 99:   L500:;
100:   }
101:   /*       RESET XADJ TO ITS CORRECT SIGN.*/
102:   i__1 = mp1end;
103:   for (i = mp1beg; i <= i__1; ++i) {
104:     node       = ls[i];
105:     xadj[node] = -xadj[node];
106:   }
107:   PetscFunctionReturn(PETSC_SUCCESS);
108: }