Actual source code: fn1wd.c
1: #define PETSCMAT_DLL
3: /* fn1wd.f -- translated by f2c (version 19931217).*/
5: #include ../src/mat/order/order.h
6: EXTERN PetscErrorCode SPARSEPACKfnroot(PetscInt*, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *, PetscInt *);
8: /*****************************************************************/
9: /******** FN1WD ..... FIND ONE-WAY DISSECTORS *********/
10: /*****************************************************************/
11: /* PURPOSE - THIS SUBROUTINE FINDS ONE-WAY DISSECTORS OF */
12: /* A CONNECTED COMPONENT SPECIFIED BY MASK AND ../../... */
13: /* */
14: /* INPUT PARAMETERS - */
15: /* ../../.. - A NODE THAT DEFINES (ALONG WITH MASK) THE */
16: /* COMPONENT TO BE PROCESSED. */
17: /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE. */
18: /* */
19: /* OUTPUT PARAMETERS - */
20: /* NSEP - NUMBER OF NODES IN THE ONE-WAY DISSECTORS. */
21: /* SEP - VECTOR CONTAINING THE DISSECTOR NODES. */
22: /* */
23: /* UPDATED PARAMETER - */
24: /* MASK - NODES IN THE DISSECTOR HAVE THEIR MASK VALUES */
25: /* SET TO ZERO. */
26: /* */
27: /* WORKING PARAMETERS- */
28: /* (XLS, LS) - LEVEL STRUCTURE USED BY THE ROUTINE FN../../... */
29: /* */
30: /* PROGRAM SUBROUTINE - */
31: /* FN../../... */
32: /*****************************************************************/
35: PetscErrorCode SPARSEPACKfn1wd(PetscInt *root, PetscInt *xadj, PetscInt *adjncy,
36: PetscInt *mask, PetscInt *nsep, PetscInt *sep, PetscInt *nlvl, PetscInt *
37: xls, PetscInt *ls)
38: {
39: /* System generated locals */
40: PetscInt i__1, i__2;
42: /* Local variables */
43: PetscInt node, i, j, k;
44: PetscReal width, fnlvl;
45: PetscInt kstop, kstrt, lp1beg, lp1end;
46: PetscReal deltp1;
47: PetscInt lvlbeg, lvlend;
48: PetscInt nbr, lvl;
51: /* Parameter adjustments */
52: --ls;
53: --xls;
54: --sep;
55: --mask;
56: --adjncy;
57: --xadj;
59: SPARSEPACKfnroot(root, &xadj[1], &adjncy[1], &mask[1], nlvl, &xls[1], &ls[1]);
60: fnlvl = (PetscReal) (*nlvl);
61: *nsep = xls[*nlvl + 1] - 1;
62: width = (PetscReal) (*nsep) / fnlvl;
63: deltp1 = sqrt((width * 3.f + 13.f) / 2.f) + 1.f;
64: if (*nsep >= 50 && deltp1 <= fnlvl * .5f) {
65: goto L300;
66: }
67: /* THE COMPONENT IS TOO SMALL, OR THE LEVEL STRUCTURE */
68: /* IS VERY LONG AND NARROW. RETURN THE WHOLE COMPONENT.*/
69: i__1 = *nsep;
70: for (i = 1; i <= i__1; ++i) {
71: node = ls[i];
72: sep[i] = node;
73: mask[node] = 0;
74: }
75: return(0);
76: /* FIND THE PARALLEL DISSECTORS.*/
77: L300:
78: *nsep = 0;
79: i = 0;
80: L400:
81: ++i;
82: lvl = (PetscInt)((PetscReal) i * deltp1 + .5f);
83: if (lvl >= *nlvl) {
84: return(0);
85: }
86: lvlbeg = xls[lvl];
87: lp1beg = xls[lvl + 1];
88: lvlend = lp1beg - 1;
89: lp1end = xls[lvl + 2] - 1;
90: i__1 = lp1end;
91: for (j = lp1beg; j <= i__1; ++j) {
92: node = ls[j];
93: xadj[node] = -xadj[node];
94: }
95: /* NODES IN LEVEL LVL ARE CHOSEN TO FORM DISSECTOR. */
96: /* INCLUDE ONLY THOSE WITH NEIGHBORS IN LVL+1 LEVEL. */
97: /* XADJ IS USED TEMPORARILY TO MARK NODES IN LVL+1. */
98: i__1 = lvlend;
99: for (j = lvlbeg; j <= i__1; ++j) {
100: node = ls[j];
101: kstrt = xadj[node];
102: kstop = (i__2 = xadj[node + 1], (PetscInt)PetscAbsInt(i__2)) - 1;
103: i__2 = kstop;
104: for (k = kstrt; k <= i__2; ++k) {
105: nbr = adjncy[k];
106: if (xadj[nbr] > 0) {
107: goto L600;
108: }
109: ++(*nsep);
110: sep[*nsep] = node;
111: mask[node] = 0;
112: goto L700;
113: L600:
114: ;
115: }
116: L700:
117: ;
118: }
119: i__1 = lp1end;
120: for (j = lp1beg; j <= i__1; ++j) {
121: node = ls[j];
122: xadj[node] = -xadj[node];
123: }
124: goto L400;
125: }