Actual source code: qmdqt.c

petsc-3.7.3 2016-08-01
Report Typos and Errors
  2: /* qmdqt.f -- translated by f2c (version 19931217).*/

  4: #include <petscsys.h>
  5: #include <petsc/private/matorderimpl.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(const PetscInt *root,const PetscInt *xadj,const PetscInt *inadjncy,PetscInt *marker, PetscInt *rchsze, PetscInt *rchset, PetscInt *nbrhd)
 30: {
 31:   PetscInt *adjncy = (PetscInt*)inadjncy; /* Used as temporary and reset within this function */
 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) goto L300;

 54: /*          PLACE REACH NODES INTO THE ADJACENT LIST OF NODE*/
 55:   i__1 = jstop;
 56:   for (j = jstrt; j <= i__1; ++j) {
 57:     ++irch;
 58:     adjncy[j] = rchset[irch];
 59:     if (irch >= *rchsze) goto L400;
 60:   }
 61: /*       LINK TO OTHER SPACE PROVIDED BY THE NBRHD SET.*/
 62: L300:
 63:   ilink = adjncy[jstop + 1];
 64:   node  = -ilink;
 65:   if (ilink < 0) goto L100;
 66:   ++inhd;
 67:   node              = nbrhd[inhd];
 68:   adjncy[jstop + 1] = -node;
 69:   goto L100;
 70: /*       ALL REACHABLE NODES HAVE BEEN SAVED.  END THE ADJ LIST.*/
 71: /*       ADD ../../.. TO THE NBR LIST OF EACH NODE IN THE REACH SET.*/
 72: L400:
 73:   adjncy[j + 1] = 0;
 74:   i__1          = *rchsze;
 75:   for (irch = 1; irch <= i__1; ++irch) {
 76:     node = rchset[irch];
 77:     if (marker[node] < 0) goto L600;

 79:     jstrt = xadj[node];
 80:     jstop = xadj[node + 1] - 1;
 81:     i__2  = jstop;
 82:     for (j = jstrt; j <= i__2; ++j) {
 83:       nabor = adjncy[j];
 84:       if (marker[nabor] >= 0) goto L500;
 85:       adjncy[j] = *root;
 86:       goto L600;
 87: L500:
 88:       ;
 89:     }
 90: L600:
 91:     ;
 92:   }
 93:   return(0);
 94: }