Actual source code: qmdrch.c

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

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

  6: /*****************************************************************/
  7: /**********     QMDRCH ..... QUOT MIN DEG REACH SET    ***********/
  8: /*****************************************************************/

 10: /*    PURPOSE - THIS SUBROUTINE DETERMINES THE REACHABLE SET OF*/
 11: /*       A NODE THROUGH A GIVEN SUBSET.  THE ADJACENCY STRUCTURE*/
 12: /*       IS ASSUMED TO BE STORED IN A QUOTIENT GRAPH FORMAT.*/

 14: /*    INPUT PARAMETERS -*/
 15: /*       ROOT - THE GIVEN NODE NOT IN THE SUBSET.*/
 16: /*       (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE PAIR.*/
 17: /*       DEG - THE DEGREE VECTOR.  DEG(I) LT 0 MEANS THE NODE*/
 18: /*              BELONGS TO THE GIVEN SUBSET.*/

 20: /*    OUTPUT PARAMETERS -*/
 21: /*       (RCHSZE, RCHSET) - THE REACHABLE SET.*/
 22: /*       (NHDSZE, NBRHD) - THE NEIGHBORHOOD SET.*/

 24: /*    UPDATED PARAMETERS -*/
 25: /*       MARKER - THE MARKER VECTOR FOR REACH AND NBRHD SETS.*/
 26: /*              GT 0 MEANS THE NODE IS IN REACH SET.*/
 27: /*              LT 0 MEANS THE NODE HAS BEEN MERGED WITH*/
 28: /*              OTHERS IN THE QUOTIENT OR IT IS IN NBRHD SET.*/
 29: /*****************************************************************/
 30: PetscErrorCode SPARSEPACKqmdrch(const PetscInt *root, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *deg, PetscInt *marker, PetscInt *rchsze, PetscInt *rchset, PetscInt *nhdsze, PetscInt *nbrhd)
 31: {
 32:   /* System generated locals */
 33:   PetscInt i__1, i__2;

 35:   /* Local variables */
 36:   PetscInt node, i, j, nabor, istop, jstop, istrt, jstrt;

 38:   /*       LOOP THROUGH THE NEIGHBORS OF ROOT IN THE*/
 39:   /*       QUOTIENT GRAPH.*/

 41:   PetscFunctionBegin;
 42:   /* Parameter adjustments */
 43:   --nbrhd;
 44:   --rchset;
 45:   --marker;
 46:   --deg;
 47:   --adjncy;
 48:   --xadj;

 50:   *nhdsze = 0;
 51:   *rchsze = 0;
 52:   istrt   = xadj[*root];
 53:   istop   = xadj[*root + 1] - 1;
 54:   if (istop < istrt) PetscFunctionReturn(PETSC_SUCCESS);
 55:   i__1 = istop;
 56:   for (i = istrt; i <= i__1; ++i) {
 57:     nabor = adjncy[i];
 58:     if (!nabor) PetscFunctionReturn(PETSC_SUCCESS);
 59:     if (marker[nabor] != 0) goto L600;
 60:     if (deg[nabor] < 0) goto L200;

 62:     /*                   INCLUDE NABOR INTO THE REACHABLE SET.*/
 63:     ++(*rchsze);
 64:     rchset[*rchsze] = nabor;
 65:     marker[nabor]   = 1;
 66:     goto L600;
 67:   /*                NABOR HAS BEEN ELIMINATED. FIND NODES*/
 68:   /*                REACHABLE FROM IT.*/
 69:   L200:
 70:     marker[nabor] = -1;
 71:     ++(*nhdsze);
 72:     nbrhd[*nhdsze] = nabor;
 73:   L300:
 74:     jstrt = xadj[nabor];
 75:     jstop = xadj[nabor + 1] - 1;
 76:     i__2  = jstop;
 77:     for (j = jstrt; j <= i__2; ++j) {
 78:       node  = adjncy[j];
 79:       nabor = -node;
 80:       if (node < 0) goto L300;
 81:       else if (!node) goto L600;
 82:       else goto L400;
 83:     L400:
 84:       if (marker[node] != 0) goto L500;
 85:       ++(*rchsze);
 86:       rchset[*rchsze] = node;
 87:       marker[node]    = 1;
 88:     L500:;
 89:     }
 90:   L600:;
 91:   }
 92:   PetscFunctionReturn(PETSC_SUCCESS);
 93: }