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