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