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