Actual source code: qmdqt.c

  1: #define PETSCMAT_DLL

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

 5:  #include petsc.h

  7: /***************************************************************/
  8: /********     QMDQT  ..... QUOT MIN DEG QUOT TRANSFORM  ********/
  9: /***************************************************************/

 11: /*    PURPOSE - THIS SUBROUTINE PERFORMS THE QUOTIENT GRAPH  */
 12: /*       TRANSFORMATION AFTER A NODE HAS BEEN ELIMINATED.*/

 14: /*    INPUT PARAMETERS -*/
 15: /*       ../../.. - THE NODE JUST ELIMINATED. IT BECOMES THE*/
 16: /*              REPRESENTATIVE OF THE NEW SUPERNODE.*/
 17: /*       (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE.*/
 18: /*       (RCHSZE, RCHSET) - THE REACHABLE SET OF ../../.. IN THE*/
 19: /*              OLD QUOTIENT GRAPH.*/
 20: /*       NBRHD - THE NEIGHBORHOOD SET WHICH WILL BE MERGED*/
 21: /*              WITH ../../.. TO FORM THE NEW SUPERNODE.*/
 22: /*       MARKER - THE MARKER VECTOR.*/

 24: /*    UPDATED PARAMETER -*/
 25: /*       ADJNCY - BECOMES THE ADJNCY OF THE QUOTIENT GRAPH.*/
 26: /***************************************************************/
 29: PetscErrorCode SPARSEPACKqmdqt(PetscInt *root, PetscInt *xadj, PetscInt *adjncy, 
 30:         PetscInt *marker, PetscInt *rchsze, PetscInt *rchset, PetscInt *nbrhd)
 31: {
 32:     /* System generated locals */
 33:     PetscInt i__1, i__2;

 35:     /* Local variables */
 36:     PetscInt inhd, irch, node, ilink, j, nabor, jstop, jstrt;

 39:     /* Parameter adjustments */
 40:     --nbrhd;
 41:     --rchset;
 42:     --marker;
 43:     --adjncy;
 44:     --xadj;

 46:     irch = 0;
 47:     inhd = 0;
 48:     node = *root;
 49: L100:
 50:     jstrt = xadj[node];
 51:     jstop = xadj[node + 1] - 2;
 52:     if (jstop < jstrt) {
 53:         goto L300;
 54:     }
 55: /*          PLACE REACH NODES INTO THE ADJACENT LIST OF NODE*/
 56:     i__1 = jstop;
 57:     for (j = jstrt; j <= i__1; ++j) {
 58:         ++irch;
 59:         adjncy[j] = rchset[irch];
 60:         if (irch >= *rchsze) {
 61:             goto L400;
 62:         }
 63:     }
 64: /*       LINK TO OTHER SPACE PROVIDED BY THE NBRHD SET.*/
 65: L300:
 66:     ilink = adjncy[jstop + 1];
 67:     node = -ilink;
 68:     if (ilink < 0) {
 69:         goto L100;
 70:     }
 71:     ++inhd;
 72:     node = nbrhd[inhd];
 73:     adjncy[jstop + 1] = -node;
 74:     goto L100;
 75: /*       ALL REACHABLE NODES HAVE BEEN SAVED.  END THE ADJ LIST.*/
 76: /*       ADD ../../.. TO THE NBR LIST OF EACH NODE IN THE REACH SET.*/
 77: L400:
 78:     adjncy[j + 1] = 0;
 79:     i__1 = *rchsze;
 80:     for (irch = 1; irch <= i__1; ++irch) {
 81:         node = rchset[irch];
 82:         if (marker[node] < 0) {
 83:             goto L600;
 84:         }
 85:         jstrt = xadj[node];
 86:         jstop = xadj[node + 1] - 1;
 87:         i__2 = jstop;
 88:         for (j = jstrt; j <= i__2; ++j) {
 89:             nabor = adjncy[j];
 90:             if (marker[nabor] >= 0) {
 91:                 goto L500;
 92:             }
 93:             adjncy[j] = *root;
 94:             goto L600;
 95: L500:
 96:             ;
 97:         }
 98: L600:
 99:         ;
100:     }
101:     return(0);
102: }