Actual source code: gennd.c
1: /* gennd.f -- translated by f2c (version 19931217).*/
3: #include <petscsys.h>
4: #include <petsc/private/matorderimpl.h>
6: PetscErrorCode SPARSEPACKrevrse(const PetscInt *n, PetscInt *perm)
7: {
8: /* System generated locals */
9: PetscInt i__1;
11: /* Local variables */
12: PetscInt swap, i, m, in;
14: PetscFunctionBegin;
15: /* Parameter adjustments */
16: --perm;
18: in = *n;
19: m = *n / 2;
20: i__1 = m;
21: for (i = 1; i <= i__1; ++i) {
22: swap = perm[i];
23: perm[i] = perm[in];
24: perm[in] = swap;
25: --in;
26: }
27: PetscFunctionReturn(PETSC_SUCCESS);
28: }
30: /*****************************************************************/
31: /********* GENND ..... GENERAL NESTED DISSECTION *********/
32: /*****************************************************************/
34: /* PURPOSE - SUBROUTINE GENND FINDS A NESTED DISSECTION*/
35: /* ORDERING FOR A GENERAL GRAPH.*/
37: /* INPUT PARAMETERS -*/
38: /* NEQNS - NUMBER OF EQUATIONS.*/
39: /* (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR.*/
41: /* OUTPUT PARAMETERS -*/
42: /* PERM - THE NESTED DISSECTION ORDERING.*/
44: /* WORKING PARAMETERS -*/
45: /* MASK - IS USED TO MASK OFF VARIABLES THAT HAVE*/
46: /* BEEN NUMBERED DURING THE ORDERNG PROCESS.*/
47: /* (XLS, LS) - THIS LEVEL STRUCTURE PAIR IS USED AS*/
48: /* TEMPORARY STORAGE BY FNROOT.*/
50: /* PROGRAM SUBROUTINES -*/
51: /* FNDSEP, REVRSE.*/
52: /*****************************************************************/
54: PetscErrorCode SPARSEPACKgennd(const PetscInt *neqns, const PetscInt *xadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *perm, PetscInt *xls, PetscInt *ls)
55: {
56: /* System generated locals */
57: PetscInt i__1;
59: /* Local variables */
60: PetscInt nsep, root, i;
61: PetscInt num;
63: PetscFunctionBegin;
64: /* Parameter adjustments */
65: --ls;
66: --xls;
67: --perm;
68: --mask;
69: --adjncy;
70: --xadj;
72: i__1 = *neqns;
73: for (i = 1; i <= i__1; ++i) mask[i] = 1;
74: num = 0;
75: i__1 = *neqns;
76: for (i = 1; i <= i__1; ++i) {
77: /* FOR EACH MASKED COMPONENT ...*/
78: L200:
79: if (!mask[i]) goto L300;
80: root = i;
81: /* FIND A SEPARATOR AND NUMBER THE NODES NEXT.*/
82: PetscCall(SPARSEPACKfndsep(&root, &xadj[1], &adjncy[1], &mask[1], &nsep, &perm[num + 1], &xls[1], &ls[1]));
83: num += nsep;
84: if (num >= *neqns) goto L400;
85: goto L200;
86: L300:;
87: }
88: /* SINCE SEPARATORS FOUND FIRST SHOULD BE ORDERED*/
89: /* LAST, ROUTINE REVRSE IS CALLED TO ADJUST THE*/
90: /* ORDERING VECTOR.*/
91: L400:
92: PetscCall(SPARSEPACKrevrse(neqns, &perm[1]));
93: PetscFunctionReturn(PETSC_SUCCESS);
94: }