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