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