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