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