Actual source code: degree.c
1: /*$Id: degree.c,v 1.21 2001/03/23 23:22:51 balay Exp $*/
2: /* degree.f -- translated by f2c (version 19931217).*/
4: #include petsc.h
5: #include src/mat/order/order.h
7: /*****************************************************************/
8: /********* DEGREE ..... DEGREE IN MASKED COMPONENT *********/
9: /*****************************************************************/
11: /* PURPOSE - THIS ROUTINE COMPUTES THE DEGREES OF THE NODES*/
12: /* IN THE CONNECTED COMPONENT SPECIFIED BY MASK AND ../../..*/
13: /* NODES FOR WHICH MASK IS ZERO ARE IGNORED.*/
15: /* INPUT PARAMETER -*/
16: /* ../../.. - IS THE INPUT NODE THAT DEFINES THE COMPONENT.*/
17: /* (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR.*/
18: /* MASK - SPECIFIES A SECTION SUBGRAPH.*/
20: /* OUTPUT PARAMETERS -*/
21: /* DEG - ARRAY CONTAINING THE DEGREES OF THE NODES IN*/
22: /* THE COMPONENT.*/
23: /* CCSIZE-SIZE OF THE COMPONENT SPECIFED BY MASK AND ../../..*/
24: /* WORKING PARAMETER -*/
25: /* LS - A TEMPORARY VECTOR USED TO STORE THE NODES OF THE*/
26: /* COMPONENT LEVEL BY LEVEL.*/
27: /*****************************************************************/
28: int SPARSEPACKdegree(int *root, int *xadj,int *adjncy,int *mask,int *deg,int *ccsize,int *ls)
29: {
30: /* System generated locals */
31: int i__1,i__2;
33: /* Local variables */
34: int ideg,node,i,j,jstop,jstrt,lbegin,lvlend,lvsize,
35: 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], (int)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] == 0) {
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: }