Actual source code: degree.c

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

  4: #include <../src/mat/order/order.h>

  6: /*****************************************************************/
  7: /*********     DEGREE ..... DEGREE IN MASKED COMPONENT   *********/
  8: /*****************************************************************/

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

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

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

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

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

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