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