Actual source code: degree.c

  1: #define PETSCMAT_DLL

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

 5:  #include petsc.h
 6:  #include src/mat/order/order.h

  8: /*****************************************************************/
  9: /*********     DEGREE ..... DEGREE IN MASKED COMPONENT   *********/
 10: /*****************************************************************/

 12: /*    PURPOSE - THIS ROUTINE COMPUTES THE DEGREES OF THE NODES*/
 13: /*       IN THE CONNECTED COMPONENT SPECIFIED BY MASK AND ../../..*/
 14: /*       NODES FOR WHICH MASK IS ZERO ARE IGNORED.*/

 16: /*    INPUT PARAMETER -*/
 17: /*       ../../.. - IS THE INPUT NODE THAT DEFINES THE COMPONENT.*/
 18: /*       (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR.*/
 19: /*       MASK - SPECIFIES A SECTION SUBGRAPH.*/

 21: /*    OUTPUT PARAMETERS -*/
 22: /*       DEG - ARRAY CONTAINING THE DEGREES OF THE NODES IN*/
 23: /*             THE COMPONENT.*/
 24: /*       CCSIZE-SIZE OF THE COMPONENT SPECIFED BY MASK AND ../../..*/
 25: /*    WORKING PARAMETER -*/
 26: /*       LS - A TEMPORARY VECTOR USED TO STORE THE NODES OF THE*/
 27: /*              COMPONENT LEVEL BY LEVEL.*/
 28: /*****************************************************************/
 31: PetscErrorCode SPARSEPACKdegree(PetscInt *root, PetscInt *xadj,PetscInt *adjncy,PetscInt *mask,PetscInt *deg,PetscInt *ccsize,PetscInt *ls)
 32: {
 33:     /* System generated locals */
 34:     PetscInt i__1,i__2;

 36:     /* Local variables */
 37:     PetscInt ideg,node,i,j,jstop,jstrt,lbegin,lvlend,lvsize,nbr;
 38: /*       INITIALIZATION ...*/
 39: /*       THE ARRAY XADJ IS USED AS A TEMPORARY MARKER TO*/
 40: /*       INDICATE WHICH NODES HAVE BEEN CONSIDERED SO FAR.*/

 43:     /* Parameter adjustments */
 44:     --ls;
 45:     --deg;
 46:     --mask;
 47:     --adjncy;
 48:     --xadj;

 50:     ls[1] = *root;
 51:     xadj[*root] = -xadj[*root];
 52:     lvlend = 0;
 53:     *ccsize = 1;
 54: /*       LBEGIN IS THE POINTER TO THE BEGINNING OF THE CURRENT*/
 55: /*       LEVEL, AND LVLEND POINTS TO THE END OF THIS LEVEL.*/
 56: L100:
 57:     lbegin = lvlend + 1;
 58:     lvlend = *ccsize;
 59: /*       FIND THE DEGREES OF NODES IN THE CURRENT LEVEL,*/
 60: /*       AND AT THE SAME TIME, GENERATE THE NEXT LEVEL.*/
 61:     i__1 = lvlend;
 62:     for (i = lbegin; i <= i__1; ++i) {
 63:         node = ls[i];
 64:         jstrt = -xadj[node];
 65:         jstop = (i__2 = xadj[node + 1], (PetscInt)PetscAbsInt(i__2)) - 1;
 66:         ideg = 0;
 67:         if (jstop < jstrt) {
 68:             goto L300;
 69:         }
 70:         i__2 = jstop;
 71:         for (j = jstrt; j <= i__2; ++j) {
 72:             nbr = adjncy[j];
 73:             if (!mask[nbr]) {
 74:                 goto L200;
 75:             }
 76:             ++ideg;
 77:             if (xadj[nbr] < 0) {
 78:                 goto L200;
 79:             }
 80:             xadj[nbr] = -xadj[nbr];
 81:             ++(*ccsize);
 82:             ls[*ccsize] = nbr;
 83: L200:
 84:             ;
 85:         }
 86: L300:
 87:         deg[node] = ideg;
 88:     }
 89: /*       COMPUTE THE CURRENT LEVEL WIDTH. */
 90: /*       IF IT IS NONZERO, GENERATE ANOTHER LEVEL.*/
 91:     lvsize = *ccsize - lvlend;
 92:     if (lvsize > 0) {
 93:         goto L100;
 94:     }
 95: /*       RESET XADJ TO ITS CORRECT SIGN AND RETURN. */
 96: /*       ------------------------------------------*/
 97:     i__1 = *ccsize;
 98:     for (i = 1; i <= i__1; ++i) {
 99:         node = ls[i];
100:         xadj[node] = -xadj[node];
101:     }
102:     return(0);
103: }