2: /* degree.f -- translated by f2c (version 19931217).*/
4: #include <petsc/private/matorderimpl.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: /*****************************************************************/
27: PetscErrorCode SPARSEPACKdegree(const PetscInt *root,const PetscInt *inxadj,const PetscInt *adjncy,PetscInt *mask,PetscInt *deg,PetscInt *ccsize,PetscInt *ls) 28: {
29: PetscInt *xadj = (PetscInt*)inxadj; /* Used as temporary and reset within this function */
30: /* System generated locals */
31: PetscInt i__1,i__2;
33: /* Local variables */
34: PetscInt ideg,node,i,j,jstop,jstrt,lbegin,lvlend,lvsize,nbr;
35: /* INITIALIZATION ...*/
36: /* THE ARRAY XADJ IS USED AS A TEMPORARY MARKER TO*/
37: /* INDICATE WHICH NODES HAVE BEEN CONSIDERED SO FAR.*/
40: /* Parameter adjustments */
41: --ls;
42: --deg;
43: --mask;
44: --adjncy;
45: --xadj;
47: ls[1] = *root;
48: xadj[*root] = -xadj[*root];
49: lvlend = 0;
50: *ccsize = 1;
51: /* LBEGIN IS THE POINTER TO THE BEGINNING OF THE CURRENT*/
52: /* LEVEL, AND LVLEND POINTS TO THE END OF THIS LEVEL.*/
53: L100: 54: lbegin = lvlend + 1;
55: lvlend = *ccsize;
56: /* FIND THE DEGREES OF NODES IN THE CURRENT LEVEL,*/
57: /* AND AT THE SAME TIME, GENERATE THE NEXT LEVEL.*/
58: i__1 = lvlend;
59: for (i = lbegin; i <= i__1; ++i) {
60: node = ls[i];
61: jstrt = -xadj[node];
62: i__2 = xadj[node + 1];
63: jstop = (PetscInt)PetscAbsInt(i__2) - 1;
64: ideg = 0;
65: if (jstop < jstrt) goto L300;
66: i__2 = jstop;
67: for (j = jstrt; j <= i__2; ++j) {
68: nbr = adjncy[j];
69: if (!mask[nbr]) goto L200;
70: ++ideg;
71: if (xadj[nbr] < 0) goto L200;
72: xadj[nbr] = -xadj[nbr];
73: ++(*ccsize);
74: ls[*ccsize] = nbr;
75: L200: 76: ;
77: }
78: L300: 79: deg[node] = ideg;
80: }
81: /* COMPUTE THE CURRENT LEVEL WIDTH. */
82: /* IF IT IS NONZERO, GENERATE ANOTHER LEVEL.*/
83: lvsize = *ccsize - lvlend;
84: if (lvsize > 0) goto L100;
85: /* RESET XADJ TO ITS CORRECT SIGN AND RETURN. */
86: /* ------------------------------------------*/
87: i__1 = *ccsize;
88: for (i = 1; i <= i__1; ++i) {
89: node = ls[i];
90: xadj[node] = -xadj[node];
91: }
92: return(0);
93: }