Actual source code: gen1wd.c
1: /* gen1wd.f -- translated by f2c (version 19931217).*/
3: #include petsc.h
5: EXTERN PetscErrorCode SPARSEPACKfn1wd(PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
6: EXTERN PetscErrorCode SPARSEPACKrevrse(PetscInt*,PetscInt*),SPARSEPACKrootls(PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
8: /*****************************************************************/
9: /*********** GEN1WD ..... GENERAL ONE-WAY DISSECTION ********/
10: /*****************************************************************/
12: /* PURPOSE - GEN1WD FINDS A ONE-WAY DISSECTION PARTITIONING*/
13: /* FOR A GENERAL GRAPH. FN1WD IS USED FOR EACH CONNECTED*/
14: /* COMPONENT.*/
16: /* INPUT PARAMETERS -*/
17: /* NEQNS - NUMBER OF EQUATIONS.*/
18: /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE PAIR.*/
20: /* OUTPUT PARAMETERS -*/
21: /* (NBLKS, XBLK) - THE PARTITIONING FOUND.*/
22: /* PERM - THE ONE-WAY DISSECTION ORDERING.*/
24: /* WORKING VECTORS -*/
25: /* MASK - IS USED TO MARK VARIABLES THAT HAVE*/
26: /* BEEN NUMBERED DURING THE ORDERING PROCESS.*/
27: /* (XLS, LS) - LEVEL STRUCTURE USED BY ../../..LS.*/
29: /* PROGRAM SUBROUTINES -*/
30: /* FN1WD, REVRSE, ../../..LS.*/
31: /****************************************************************/
34: PetscErrorCode SPARSEPACKgen1wd(PetscInt *neqns, PetscInt *xadj, PetscInt *adjncy,
35: PetscInt *mask, PetscInt *nblks, PetscInt *xblk, PetscInt *perm, PetscInt *xls, PetscInt *ls)
36: {
37: /* System generated locals */
38: PetscInt i__1, i__2, i__3;
40: /* Local variables */
41: PetscInt node, nsep, lnum, nlvl, root;
42: PetscInt i, j, k, ccsize;
43: PetscInt num;
46: /* Parameter adjustments */
47: --ls;
48: --xls;
49: --perm;
50: --xblk;
51: --mask;
52: --xadj;
53: --adjncy;
55: i__1 = *neqns;
56: for (i = 1; i <= i__1; ++i) {
57: mask[i] = 1;
58: }
59: *nblks = 0;
60: num = 0;
61: i__1 = *neqns;
62: for (i = 1; i <= i__1; ++i) {
63: if (!mask[i]) {
64: goto L400;
65: }
66: /* FIND A ONE-WAY DISSECTOR FOR EACH COMPONENT.*/
67: root = i;
68: SPARSEPACKfn1wd(&root, &xadj[1], &adjncy[1], &mask[1], &nsep, &perm[num + 1], &
69: nlvl, &xls[1], &ls[1]);
70: num += nsep;
71: ++(*nblks);
72: xblk[*nblks] = *neqns - num + 1;
73: ccsize = xls[nlvl + 1] - 1;
74: /* NUMBER THE REMAINING NODES IN THE COMPONENT.*/
75: /* EACH COMPONENT IN THE REMAINING SUBGRAPH FORMS*/
76: /* A NEW BLOCK IN THE PARTITIONING.*/
77: i__2 = ccsize;
78: for (j = 1; j <= i__2; ++j) {
79: node = ls[j];
80: if (!mask[node]) {
81: goto L300;
82: }
83: SPARSEPACKrootls(&node, &xadj[1], &adjncy[1], &mask[1], &nlvl, &xls[1], &
84: perm[num + 1]);
85: lnum = num + 1;
86: num = num + xls[nlvl + 1] - 1;
87: ++(*nblks);
88: xblk[*nblks] = *neqns - num + 1;
89: i__3 = num;
90: for (k = lnum; k <= i__3; ++k) {
91: node = perm[k];
92: mask[node] = 0;
93: }
94: if (num > *neqns) {
95: goto L500;
96: }
97: L300:
98: ;
99: }
100: L400:
101: ;
102: }
103: /* SINCE DISSECTORS FOUND FIRST SHOULD BE ORDERED LAST,*/
104: /* ROUTINE REVRSE IS CALLED TO ADJUST THE ORDERING*/
105: /* VECTOR, AND THE BLOCK INDEX VECTOR.*/
106: L500:
107: SPARSEPACKrevrse(neqns, &perm[1]);
108: SPARSEPACKrevrse(nblks, &xblk[1]);
109: xblk[*nblks + 1] = *neqns + 1;
110: return(0);
111: }