Actual source code: rootls.c
1: #define PETSCMAT_DLL
3: /* rootls.f -- translated by f2c (version 19931217).*/
5: #include petscsys.h
7: /*****************************************************************/
8: /********* ../../..LS ..... ../../..ED LEVEL STRUCTURE **********/
9: /*****************************************************************/
10: /* PURPOSE - ../../..LS GENERATES THE LEVEL STRUCTURE ../../..ED */
11: /* AT THE INPUT NODE CALLED ../../... ONLY THOSE NODES FOR*/
12: /* WHICH MASK IS NONZERO WILL BE CONSIDERED.*/
13: /* */
14: /* INPUT PARAMETERS - */
15: /* ../../.. - THE NODE AT WHICH THE LEVEL STRUCTURE IS TO*/
16: /* BE ../../..ED.*/
17: /* (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR FOR THE*/
18: /* GIVEN GRAPH.*/
19: /* MASK - IS USED TO SPECIFY A SECTION SUBGRAPH. NODES*/
20: /* WITH MASK(I)=0 ARE IGNORED.*/
21: /* OUTPUT PARAMETERS -*/
22: /* NLVL - IS THE NUMBER OF LEVELS IN THE LEVEL STRUCTURE.*/
23: /* (XLS, LS) - ARRAY PAIR FOR THE ../../..ED LEVEL STRUCTURE.*/
24: /*****************************************************************/
27: PetscErrorCode SPARSEPACKrootls(PetscInt *root, PetscInt *xadj, PetscInt *adjncy,
28: PetscInt *mask, PetscInt *nlvl, PetscInt *xls, PetscInt *ls)
29: {
30: /* System generated locals */
31: PetscInt i__1, i__2;
33: /* Local variables */
34: PetscInt node, i, j, jstop, jstrt, lbegin, ccsize, lvlend, lvsize,
35: nbr;
37: /* INITIALIZATION ...*/
41: /* Parameter adjustments */
42: --ls;
43: --xls;
44: --mask;
45: --adjncy;
46: --xadj;
48: mask[*root] = 0;
49: ls[1] = *root;
50: *nlvl = 0;
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: L200:
56: lbegin = lvlend + 1;
57: lvlend = ccsize;
58: ++(*nlvl);
59: xls[*nlvl] = lbegin;
60: /* GENERATE THE NEXT LEVEL BY FINDING ALL THE MASKED */
61: /* NEIGHBORS OF NODES IN THE CURRENT LEVEL.*/
62: i__1 = lvlend;
63: for (i = lbegin; i <= i__1; ++i) {
64: node = ls[i];
65: jstrt = xadj[node];
66: jstop = xadj[node + 1] - 1;
67: if (jstop < jstrt) {
68: goto L400;
69: }
70: i__2 = jstop;
71: for (j = jstrt; j <= i__2; ++j) {
72: nbr = adjncy[j];
73: if (!mask[nbr]) {
74: goto L300;
75: }
76: ++ccsize;
77: ls[ccsize] = nbr;
78: mask[nbr] = 0;
79: L300:
80: ;
81: }
82: L400:
83: ;
84: }
85: /* COMPUTE THE CURRENT LEVEL WIDTH.*/
86: /* IF IT IS NONZERO, GENERATE THE NEXT LEVEL.*/
87: lvsize = ccsize - lvlend;
88: if (lvsize > 0) {
89: goto L200;
90: }
91: /* RESET MASK TO ONE FOR THE NODES IN THE LEVEL STRUCTURE.*/
92: xls[*nlvl + 1] = lvlend + 1;
93: i__1 = ccsize;
94: for (i = 1; i <= i__1; ++i) {
95: node = ls[i];
96: mask[node] = 1;
97: }
98: return(0);
99: }