Actual source code: gennd.c


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

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

  7: PetscErrorCode SPARSEPACKrevrse(const PetscInt *n,PetscInt *perm)
  8: {
  9:   /* System generated locals */
 10:   PetscInt i__1;

 12:   /* Local variables */
 13:   PetscInt swap,i,m,in;

 16:   /* Parameter adjustments */
 17:   --perm;

 19:   in   = *n;
 20:   m    = *n / 2;
 21:   i__1 = m;
 22:   for (i = 1; i <= i__1; ++i) {
 23:     swap     = perm[i];
 24:     perm[i]  = perm[in];
 25:     perm[in] = swap;
 26:     --in;
 27:   }
 28:   return(0);
 29: }

 31: /*****************************************************************/
 32: /*********     GENND ..... GENERAL NESTED DISSECTION     *********/
 33: /*****************************************************************/

 35: /*    PURPOSE - SUBROUTINE GENND FINDS A NESTED DISSECTION*/
 36: /*       ORDERING FOR A GENERAL GRAPH.*/

 38: /*    INPUT PARAMETERS -*/
 39: /*       NEQNS - NUMBER OF EQUATIONS.*/
 40: /*       (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR.*/

 42: /*    OUTPUT PARAMETERS -*/
 43: /*       PERM - THE NESTED DISSECTION ORDERING.*/

 45: /*    WORKING PARAMETERS -*/
 46: /*       MASK - IS USED TO MASK OFF VARIABLES THAT HAVE*/
 47: /*              BEEN NUMBERED DURING THE ORDERNG PROCESS.*/
 48: /*       (XLS, LS) - THIS LEVEL STRUCTURE PAIR IS USED AS*/
 49: /*              TEMPORARY STORAGE BY FN../../...*/

 51: /*    PROGRAM SUBROUTINES -*/
 52: /*       FNDSEP, REVRSE.*/
 53: /*****************************************************************/

 55: PetscErrorCode SPARSEPACKgennd(const PetscInt *neqns,const PetscInt *xadj,const PetscInt *adjncy,PetscInt *mask,PetscInt *perm,PetscInt *xls,PetscInt *ls)
 56: {
 57:   /* System generated locals */
 58:   PetscInt i__1;

 60:   /* Local variables */
 61:   PetscInt nsep,root,i;
 62:   PetscInt num;

 65:   /* Parameter adjustments */
 66:   --ls;
 67:   --xls;
 68:   --perm;
 69:   --mask;
 70:   --adjncy;
 71:   --xadj;

 73:   i__1 = *neqns;
 74:   for (i = 1; i <= i__1; ++i) mask[i] = 1;
 75:   num  = 0;
 76:   i__1 = *neqns;
 77:   for (i = 1; i <= i__1; ++i) {
 78: /*           FOR EACH MASKED COMPONENT ...*/
 79: L200:
 80:     if (!mask[i]) goto L300;
 81:     root = i;
 82: /*              FIND A SEPARATOR AND NUMBER THE NODES NEXT.*/
 83:     SPARSEPACKfndsep(&root,&xadj[1],&adjncy[1],&mask[1],&nsep,&perm[num + 1],&xls[1],&ls[1]);
 84:     num += nsep;
 85:     if (num >= *neqns) goto L400;
 86:     goto L200;
 87: L300:
 88:     ;
 89:   }
 90: /*        SINCE SEPARATORS FOUND FIRST SHOULD BE ORDERED*/
 91: /*        LAST, ROUTINE REVRSE IS CALLED TO ADJUST THE*/
 92: /*        ORDERING VECTOR.*/
 93: L400:
 94:   SPARSEPACKrevrse(neqns,&perm[1]);
 95:   return(0);
 96: }