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;

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

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

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

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

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

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

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

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

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

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

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

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