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