Actual source code: genqmd.c


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

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

  7: /******************************************************************/
  8: /***********    GENQMD ..... QUOT MIN DEGREE ORDERING    **********/
  9: /******************************************************************/
 10: /*    PURPOSE - THIS ROUTINE IMPLEMENTS THE MINIMUM DEGREE        */
 11: /*       ALGORITHM.  IT MAKES USE OF THE IMPLICIT REPRESENT-      */
 12: /*       ATION OF THE ELIMINATION GRAPHS BY QUOTIENT GRAPHS,      */
 13: /*       AND THE NOTION OF INDISTINGUISHABLE NODES.               */
 14: /*       CAUTION - THE ADJACENCY VECTOR ADJNCY WILL BE            */
 15: /*       DESTROYED.                                               */
 16: /*                                                                */
 17: /*    INPUT PARAMETERS -                                          */
 18: /*       NEQNS - NUMBER OF EQUATIONS.                             */
 19: /*       (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE.                */
 20: /*                                                                */
 21: /*    OUTPUT PARAMETERS -                                         */
 22: /*       PERM - THE MINIMUM DEGREE ORDERING.                      */
 23: /*       INVP - THE INVERSE OF PERM.                              */
 24: /*                                                                */
 25: /*    WORKING PARAMETERS -                                        */
 26: /*       DEG - THE DEGREE VECTOR. DEG(I) IS NEGATIVE MEANS        */
 27: /*              NODE I HAS BEEN NUMBERED.                         */
 28: /*       MARKER - A MARKER VECTOR, WHERE MARKER(I) IS             */
 29: /*              NEGATIVE MEANS NODE I HAS BEEN MERGED WITH        */
 30: /*              ANOTHER NODE AND THUS CAN BE IGNORED.             */
 31: /*       RCHSET - VECTOR USED FOR THE REACHABLE SET.              */
 32: /*       NBRHD - VECTOR USED FOR THE NEIGHBORHOOD SET.            */
 33: /*       QSIZE - VECTOR USED TO STORE THE SIZE OF                 */
 34: /*              INDISTINGUISHABLE SUPERNODES.                     */
 35: /*       QLINK - VECTOR TO STORE INDISTINGUISHABLE NODES,         */
 36: /*              I, QLINK(I), QLINK(QLINK(I)) ... ARE THE          */
 37: /*              MEMBERS OF THE SUPERNODE REPRESENTED BY I.        */
 38: /*                                                                */
 39: /*    PROGRAM SUBROUTINES -                                       */
 40: /*       QMDRCH, QMDQT, QMDUPD.                                   */
 41: /*                                                                */
 42: /******************************************************************/
 43: /*                                                                */
 44: /*                                                                */
 45: PetscErrorCode SPARSEPACKgenqmd(const PetscInt *neqns, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *perm, PetscInt *invp, PetscInt *deg, PetscInt *marker, PetscInt *rchset, PetscInt *nbrhd, PetscInt *qsize, PetscInt *qlink, PetscInt *nofsub)
 46: {
 47:   /* System generated locals */
 48:   PetscInt i__1;

 50:   /* Local variables */
 51:   PetscInt ndeg, irch, node, nump1, j, inode;
 52:   PetscInt ip, np, mindeg, search;
 53:   PetscInt nhdsze, nxnode, rchsze, thresh, num;

 55:   /*       INITIALIZE DEGREE VECTOR AND OTHER WORKING VARIABLES.   */

 57:   PetscFunctionBegin;
 58:   /* Parameter adjustments */
 59:   --qlink;
 60:   --qsize;
 61:   --nbrhd;
 62:   --rchset;
 63:   --marker;
 64:   --deg;
 65:   --invp;
 66:   --perm;
 67:   --adjncy;
 68:   --xadj;

 70:   mindeg  = *neqns;
 71:   *nofsub = 0;
 72:   i__1    = *neqns;
 73:   for (node = 1; node <= i__1; ++node) {
 74:     perm[node]   = node;
 75:     invp[node]   = node;
 76:     marker[node] = 0;
 77:     qsize[node]  = 1;
 78:     qlink[node]  = 0;
 79:     ndeg         = xadj[node + 1] - xadj[node];
 80:     deg[node]    = ndeg;
 81:     if (ndeg < mindeg) mindeg = ndeg;
 82:   }
 83:   num = 0;
 84: /*       PERFORM THRESHOLD SEARCH TO GET A NODE OF MIN DEGREE.   */
 85: /*       VARIABLE SEARCH POINTS TO WHERE SEARCH SHOULD START.    */
 86: L200:
 87:   search = 1;
 88:   thresh = mindeg;
 89:   mindeg = *neqns;
 90: L300:
 91:   nump1 = num + 1;
 92:   if (nump1 > search) search = nump1;
 93:   i__1 = *neqns;
 94:   for (j = search; j <= i__1; ++j) {
 95:     node = perm[j];
 96:     if (marker[node] < 0) goto L400;
 97:     ndeg = deg[node];
 98:     if (ndeg <= thresh) goto L500;
 99:     if (ndeg < mindeg) mindeg = ndeg;
100:   L400:;
101:   }
102:   goto L200;
103: /*          NODE HAS MINIMUM DEGREE. FIND ITS REACHABLE SETS BY    */
104: /*          CALLING QMDRCH.                                        */
105: L500:
106:   search = j;
107:   *nofsub += deg[node];
108:   marker[node] = 1;
109:   PetscCall(SPARSEPACKqmdrch(&node, &xadj[1], &adjncy[1], &deg[1], &marker[1], &rchsze, &rchset[1], &nhdsze, &nbrhd[1]));
110:   /*          ELIMINATE ALL NODES INDISTINGUISHABLE FROM NODE.       */
111:   /*          THEY ARE GIVEN BY NODE, QLINK(NODE), ....              */
112:   nxnode = node;
113: L600:
114:   ++num;
115:   np           = invp[nxnode];
116:   ip           = perm[num];
117:   perm[np]     = ip;
118:   invp[ip]     = np;
119:   perm[num]    = nxnode;
120:   invp[nxnode] = num;
121:   deg[nxnode]  = -1;
122:   nxnode       = qlink[nxnode];
123:   if (nxnode > 0) goto L600;
124:   if (rchsze <= 0) goto L800;

126:   /*             UPDATE THE DEGREES OF THE NODES IN THE REACHABLE     */
127:   /*             SET AND IDENTIFY INDISTINGUISHABLE NODES.            */
128:   PetscCall(SPARSEPACKqmdupd(&xadj[1], &adjncy[1], &rchsze, &rchset[1], &deg[1], &qsize[1], &qlink[1], &marker[1], &rchset[rchsze + 1], &nbrhd[nhdsze + 1]));

130:   /*             RESET MARKER VALUE OF NODES IN REACH SET.            */
131:   /*             UPDATE THRESHOLD VALUE FOR CYCLIC SEARCH.            */
132:   /*             ALSO CALL QMDQT TO FORM NEW QUOTIENT GRAPH.          */
133:   marker[node] = 0;
134:   i__1         = rchsze;
135:   for (irch = 1; irch <= i__1; ++irch) {
136:     inode = rchset[irch];
137:     if (marker[inode] < 0) goto L700;

139:     marker[inode] = 0;
140:     ndeg          = deg[inode];
141:     if (ndeg < mindeg) mindeg = ndeg;
142:     if (ndeg > thresh) goto L700;
143:     mindeg = thresh;
144:     thresh = ndeg;
145:     search = invp[inode];
146:   L700:;
147:   }
148:   if (nhdsze > 0) PetscCall(SPARSEPACKqmdqt(&node, &xadj[1], &adjncy[1], &marker[1], &rchsze, &rchset[1], &nbrhd[1]));
149: L800:
150:   if (num < *neqns) goto L300;
151:   PetscFunctionReturn(PETSC_SUCCESS);
152: }