Actual source code: qmdrch.c

  1: #define PETSCMAT_DLL

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

 5:  #include petsc.h

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

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

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

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

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

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

 43: /*       LOOP THROUGH THE NEIGHBORS OF ../../.. IN THE*/
 44: /*       QUOTIENT GRAPH.*/


 48:     /* Parameter adjustments */
 49:     --nbrhd;
 50:     --rchset;
 51:     --marker;
 52:     --deg;
 53:     --adjncy;
 54:     --xadj;

 56:     *nhdsze = 0;
 57:     *rchsze = 0;
 58:     istrt = xadj[*root];
 59:     istop = xadj[*root + 1] - 1;
 60:     if (istop < istrt) {
 61:         return(0);
 62:     }
 63:     i__1 = istop;
 64:     for (i = istrt; i <= i__1; ++i) {
 65:         nabor = adjncy[i];
 66:         if (!nabor) {
 67:             return(0);
 68:         }
 69:         if (marker[nabor] != 0) {
 70:             goto L600;
 71:         }
 72:         if (deg[nabor] < 0) {
 73:             goto L200;
 74:         }
 75: /*                   INCLUDE NABOR INTO THE REACHABLE SET.*/
 76:         ++(*rchsze);
 77:         rchset[*rchsze] = nabor;
 78:         marker[nabor] = 1;
 79:         goto L600;
 80: /*                NABOR HAS BEEN ELIMINATED. FIND NODES*/
 81: /*                REACHABLE FROM IT.*/
 82: L200:
 83:         marker[nabor] = -1;
 84:         ++(*nhdsze);
 85:         nbrhd[*nhdsze] = nabor;
 86: L300:
 87:         jstrt = xadj[nabor];
 88:         jstop = xadj[nabor + 1] - 1;
 89:         i__2 = jstop;
 90:         for (j = jstrt; j <= i__2; ++j) {
 91:             node = adjncy[j];
 92:             nabor = -node;
 93:             if (node < 0) {
 94:                 goto L300;
 95:             } else if (!node) {
 96:                 goto L600;
 97:             } else {
 98:                 goto L400;
 99:             }
100: L400:
101:             if (marker[node] != 0) {
102:                 goto L500;
103:             }
104:             ++(*rchsze);
105:             rchset[*rchsze] = node;
106:             marker[node] = 1;
107: L500:
108:             ;
109:         }
110: L600:
111:         ;
112:     }
113:     return(0);
114: }