Actual source code: fndsep.c
petsc-3.9.4 2018-09-11
2: /* fndsep.f -- translated by f2c (version 19931217).
3: */
5: #include <petsc/private/matorderimpl.h>
7: /*****************************************************************/
8: /************* FNDSEP ..... FIND SEPARATOR *************/
9: /*****************************************************************/
10: /* PURPOSE - THIS ROUTINE IS USED TO FIND A SMALL */
11: /* SEPARATOR FOR A CONNECTED COMPONENT SPECIFIED */
12: /* BY MASK IN THE GIVEN GRAPH. */
13: /* */
14: /* INPUT PARAMETERS - */
15: /* ../../.. - IS THE NODE THAT DETERMINES THE MASKED */
16: /* COMPONENT. */
17: /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE PAIR. */
18: /* */
19: /* OUTPUT PARAMETERS - */
20: /* NSEP - NUMBER OF VARIABLES IN THE SEPARATOR. */
21: /* SEP - VECTOR CONTAINING THE SEPARATOR NODES. */
22: /* */
23: /* UPDATED PARAMETER - */
24: /* MASK - NODES IN THE SEPARATOR HAVE THEIR MASK */
25: /* VALUES SET TO ZERO. */
26: /* */
27: /* WORKING PARAMETERS - */
28: /* (XLS, LS) - LEVEL STRUCTURE PAIR FOR LEVEL STRUCTURE */
29: /* FOUND BY FN../../... */
30: /* */
31: /* PROGRAM SUBROUTINES - */
32: /* FN../../... */
33: /* */
34: /*****************************************************************/
35: PetscErrorCode SPARSEPACKfndsep(PetscInt *root,const PetscInt *inxadj,const PetscInt *adjncy,PetscInt *mask, PetscInt *nsep, PetscInt *sep, PetscInt *xls, PetscInt *ls)
36: {
37: PetscInt *xadj = (PetscInt*)inxadj; /* Used as temporary and reset within this function */
39: /* System generated locals */
40: PetscInt i__1, i__2;
42: /* Local variables */
43: PetscInt node, nlvl, i, j, jstop, jstrt, mp1beg, mp1end, midbeg, midend, midlvl;
44: PetscInt nbr;
47: /* Parameter adjustments */
48: --ls;
49: --xls;
50: --sep;
51: --mask;
52: --adjncy;
53: --xadj;
55: SPARSEPACKfnroot(root, &xadj[1], &adjncy[1], &mask[1], &nlvl, &xls[1], &ls[1]);
56: /* IF THE NUMBER OF LEVELS IS LESS THAN 3, RETURN */
57: /* THE WHOLE COMPONENT AS THE SEPARATOR.*/
58: if (nlvl >= 3) goto L200;
59: *nsep = xls[nlvl + 1] - 1;
60: i__1 = *nsep;
61: for (i = 1; i <= i__1; ++i) {
62: node = ls[i];
63: sep[i] = node;
64: mask[node] = 0;
65: }
66: return(0);
67: /* FIND THE MIDDLE LEVEL OF THE ../../..ED LEVEL STRUCTURE.*/
68: L200:
69: midlvl = (nlvl + 2) / 2;
70: midbeg = xls[midlvl];
71: mp1beg = xls[midlvl + 1];
72: midend = mp1beg - 1;
73: mp1end = xls[midlvl + 2] - 1;
74: /* THE SEPARATOR IS OBTAINED BY INCLUDING ONLY THOSE*/
75: /* MIDDLE-LEVEL NODES WITH NEIGHBORS IN THE MIDDLE+1*/
76: /* LEVEL. XADJ IS USED TEMPORARILY TO MARK THOSE*/
77: /* NODES IN THE MIDDLE+1 LEVEL.*/
78: i__1 = mp1end;
79: for (i = mp1beg; i <= i__1; ++i) {
80: node = ls[i];
81: xadj[node] = -xadj[node];
82: }
83: *nsep = 0;
84: i__1 = midend;
85: for (i = midbeg; i <= i__1; ++i) {
86: node = ls[i];
87: jstrt = xadj[node];
88: i__2 = xadj[node + 1];
89: jstop = (PetscInt)PetscAbsInt(i__2) - 1;
90: i__2 = jstop;
91: for (j = jstrt; j <= i__2; ++j) {
92: nbr = adjncy[j];
93: if (xadj[nbr] > 0) goto L400;
94: ++(*nsep);
95: sep[*nsep] = node;
96: mask[node] = 0;
97: goto L500;
98: L400:
99: ;
100: }
101: L500:
102: ;
103: }
104: /* RESET XADJ TO ITS CORRECT SIGN.*/
105: i__1 = mp1end;
106: for (i = mp1beg; i <= i__1; ++i) {
107: node = ls[i];
108: xadj[node] = -xadj[node];
109: }
110: return(0);
111: }