Actual source code: gennd.c
1: #define PETSCMAT_DLL
3: /* gennd.f -- translated by f2c (version 19931217).*/
5: #include petsc.h
9: PetscErrorCode SPARSEPACKrevrse(PetscInt *n,PetscInt *perm)
10: {
11: /* System generated locals */
12: PetscInt i__1;
14: /* Local variables */
15: PetscInt swap,i,m,in;
18: /* Parameter adjustments */
19: --perm;
21: in = *n;
22: m = *n / 2;
23: i__1 = m;
24: for (i = 1; i <= i__1; ++i) {
25: swap = perm[i];
26: perm[i] = perm[in];
27: perm[in] = swap;
28: --in;
29: }
30: return(0);
31: }
34: /*****************************************************************/
35: /********* GENND ..... GENERAL NESTED DISSECTION *********/
36: /*****************************************************************/
38: /* PURPOSE - SUBROUTINE GENND FINDS A NESTED DISSECTION*/
39: /* ORDERING FOR A GENERAL GRAPH.*/
41: /* INPUT PARAMETERS -*/
42: /* NEQNS - NUMBER OF EQUATIONS.*/
43: /* (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR.*/
45: /* OUTPUT PARAMETERS -*/
46: /* PERM - THE NESTED DISSECTION ORDERING.*/
48: /* WORKING PARAMETERS -*/
49: /* MASK - IS USED TO MASK OFF VARIABLES THAT HAVE*/
50: /* BEEN NUMBERED DURING THE ORDERNG PROCESS.*/
51: /* (XLS, LS) - THIS LEVEL STRUCTURE PAIR IS USED AS*/
52: /* TEMPORARY STORAGE BY FN../../...*/
54: /* PROGRAM SUBROUTINES -*/
55: /* FNDSEP, REVRSE.*/
56: /*****************************************************************/
60: PetscErrorCode SPARSEPACKgennd(PetscInt *neqns,PetscInt *xadj,PetscInt *adjncy,PetscInt *mask,PetscInt *perm,PetscInt *xls,PetscInt *ls)
61: {
62: /* System generated locals */
63: PetscInt i__1;
65: /* Local variables */
66: PetscInt nsep,root,i;
67: EXTERN PetscErrorCode SPARSEPACKfndsep(PetscInt*,PetscInt *,PetscInt *,PetscInt *,PetscInt *,PetscInt *,PetscInt *,PetscInt *);
68: PetscInt num;
71: /* Parameter adjustments */
72: --ls;
73: --xls;
74: --perm;
75: --mask;
76: --adjncy;
77: --xadj;
79: i__1 = *neqns;
80: for (i = 1; i <= i__1; ++i) {
81: mask[i] = 1;
82: }
83: num = 0;
84: i__1 = *neqns;
85: for (i = 1; i <= i__1; ++i) {
86: /* FOR EACH MASKED COMPONENT ...*/
87: L200:
88: if (!mask[i]) {
89: goto L300;
90: }
91: root = i;
92: /* FIND A SEPARATOR AND NUMBER THE NODES NEXT.*/
93: SPARSEPACKfndsep(&root,&xadj[1],&adjncy[1],&mask[1],&nsep,&perm[num + 1],
94: &xls[1],&ls[1]);
95: num += nsep;
96: if (num >= *neqns) {
97: goto L400;
98: }
99: goto L200;
100: L300:
101: ;
102: }
103: /* SINCE SEPARATORS FOUND FIRST SHOULD BE ORDERED*/
104: /* LAST, ROUTINE REVRSE IS CALLED TO ADJUST THE*/
105: /* ORDERING VECTOR.*/
106: L400:
107: SPARSEPACKrevrse(neqns,&perm[1]);
108: return(0);
109: }