Actual source code: rootls.c
1: /* rootls.f -- translated by f2c (version 19931217).*/
3: #include <petscsys.h>
4: #include <petsc/private/matorderimpl.h>
6: /*****************************************************************/
7: /********* ROOTLS ..... ROOTED LEVEL STRUCTURE **********/
8: /*****************************************************************/
9: /* PURPOSE - ROOTLS GENERATES THE LEVEL STRUCTURE ROOTED */
10: /* AT THE INPUT NODE CALLED ROOT. ONLY THOSE NODES FOR*/
11: /* WHICH MASK IS NONZERO WILL BE CONSIDERED.*/
12: /* */
13: /* INPUT PARAMETERS - */
14: /* ROOT - THE NODE AT WHICH THE LEVEL STRUCTURE IS TO*/
15: /* BE ROOTED.*/
16: /* (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR FOR THE*/
17: /* GIVEN GRAPH.*/
18: /* MASK - IS USED TO SPECIFY A SECTION SUBGRAPH. NODES*/
19: /* WITH MASK(I)=0 ARE IGNORED.*/
20: /* OUTPUT PARAMETERS -*/
21: /* NLVL - IS THE NUMBER OF LEVELS IN THE LEVEL STRUCTURE.*/
22: /* (XLS, LS) - ARRAY PAIR FOR THE ROOTED LEVEL STRUCTURE.*/
23: /*****************************************************************/
24: PetscErrorCode SPARSEPACKrootls(const PetscInt *root, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nlvl, PetscInt *xls, PetscInt *ls)
25: {
26: /* System generated locals */
27: PetscInt i__1, i__2;
29: /* Local variables */
30: PetscInt node, i, j, jstop, jstrt, lbegin, ccsize, lvlend, lvsize, nbr;
32: /* INITIALIZATION ...*/
34: PetscFunctionBegin;
35: /* Parameter adjustments */
36: --ls;
37: --xls;
38: --mask;
39: --adjncy;
40: --xadj;
42: mask[*root] = 0;
43: ls[1] = *root;
44: *nlvl = 0;
45: lvlend = 0;
46: ccsize = 1;
47: /* LBEGIN IS THE POINTER TO THE BEGINNING OF THE CURRENT*/
48: /* LEVEL, AND LVLEND POINTS TO THE END OF THIS LEVEL.*/
49: L200:
50: lbegin = lvlend + 1;
51: lvlend = ccsize;
52: ++(*nlvl);
53: xls[*nlvl] = lbegin;
54: /* GENERATE THE NEXT LEVEL BY FINDING ALL THE MASKED */
55: /* NEIGHBORS OF NODES IN THE CURRENT LEVEL.*/
56: i__1 = lvlend;
57: for (i = lbegin; i <= i__1; ++i) {
58: node = ls[i];
59: jstrt = xadj[node];
60: jstop = xadj[node + 1] - 1;
61: if (jstop < jstrt) goto L400;
62: i__2 = jstop;
63: for (j = jstrt; j <= i__2; ++j) {
64: nbr = adjncy[j];
65: if (!mask[nbr]) goto L300;
66: ++ccsize;
67: ls[ccsize] = nbr;
68: mask[nbr] = 0;
69: L300:;
70: }
71: L400:;
72: }
73: /* COMPUTE THE CURRENT LEVEL WIDTH.*/
74: /* IF IT IS NONZERO, GENERATE THE NEXT LEVEL.*/
75: lvsize = ccsize - lvlend;
76: if (lvsize > 0) goto L200;
77: /* RESET MASK TO ONE FOR THE NODES IN THE LEVEL STRUCTURE.*/
78: xls[*nlvl + 1] = lvlend + 1;
79: i__1 = ccsize;
80: for (i = 1; i <= i__1; ++i) {
81: node = ls[i];
82: mask[node] = 1;
83: }
84: PetscFunctionReturn(PETSC_SUCCESS);
85: }