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: /****************************************************************/
30: int SPARSEPACKgen1wd(int *neqns, int *xadj, int *adjncy,
31: int *mask, int *nblks, int *xblk, int *perm, int *
32: xls, int *ls)
33: {
34: /* System generated locals */
35: int i__1, i__2, i__3;
37: /* Local variables */
38: int node, nsep, lnum, nlvl, root;
39: EXTERN int SPARSEPACKfn1wd(int *, int *, int *,
40: int *, int *, int *, int *, int *, int *);
41: int i, j, k, ccsize;
42: EXTERN int SPARSEPACKrevrse(int *, int *), SPARSEPACKrootls(
43: int *, int *, int *, int *, int *, int *, int *);
44: int num;
47: /* Parameter adjustments */
48: --ls;
49: --xls;
50: --perm;
51: --xblk;
52: --mask;
53: --xadj;
54: --adjncy;
56: i__1 = *neqns;
57: for (i = 1; i <= i__1; ++i) {
58: mask[i] = 1;
59: }
60: *nblks = 0;
61: num = 0;
62: i__1 = *neqns;
63: for (i = 1; i <= i__1; ++i) {
64: if (mask[i] == 0) {
65: goto L400;
66: }
67: /* FIND A ONE-WAY DISSECTOR FOR EACH COMPONENT.*/
68: root = i;
69: SPARSEPACKfn1wd(&root, &xadj[1], &adjncy[1], &mask[1], &nsep, &perm[num + 1], &
70: nlvl, &xls[1], &ls[1]);
71: num += nsep;
72: ++(*nblks);
73: xblk[*nblks] = *neqns - num + 1;
74: ccsize = xls[nlvl + 1] - 1;
75: /* NUMBER THE REMAINING NODES IN THE COMPONENT.*/
76: /* EACH COMPONENT IN THE REMAINING SUBGRAPH FORMS*/
77: /* A NEW BLOCK IN THE PARTITIONING.*/
78: i__2 = ccsize;
79: for (j = 1; j <= i__2; ++j) {
80: node = ls[j];
81: if (mask[node] == 0) {
82: goto L300;
83: }
84: SPARSEPACKrootls(&node, &xadj[1], &adjncy[1], &mask[1], &nlvl, &xls[1], &
85: perm[num + 1]);
86: lnum = num + 1;
87: num = num + xls[nlvl + 1] - 1;
88: ++(*nblks);
89: xblk[*nblks] = *neqns - num + 1;
90: i__3 = num;
91: for (k = lnum; k <= i__3; ++k) {
92: node = perm[k];
93: mask[node] = 0;
94: }
95: if (num > *neqns) {
96: goto L500;
97: }
98: L300:
99: ;
100: }
101: L400:
102: ;
103: }
104: /* SINCE DISSECTORS FOUND FIRST SHOULD BE ORDERED LAST,*/
105: /* ROUTINE REVRSE IS CALLED TO ADJUST THE ORDERING*/
106: /* VECTOR, AND THE BLOCK INDEX VECTOR.*/
107: L500:
108: SPARSEPACKrevrse(neqns, &perm[1]);
109: SPARSEPACKrevrse(nblks, &xblk[1]);
110: xblk[*nblks + 1] = *neqns + 1;
111: return(0);
112: }