Actual source code: genrcm.c
1: /* genrcm.f -- translated by f2c (version 19931217).*/
3: #include petsc.h
5: /*****************************************************************/
6: /*****************************************************************/
7: /********* GENRCM ..... GENERAL REVERSE CUTHILL MCKEE ********/
8: /*****************************************************************/
10: /* PURPOSE - GENRCM FINDS THE REVERSE CUTHILL-MCKEE*/
11: /* ORDERING FOR A GENERAL GRAPH. FOR EACH CONNECTED*/
12: /* COMPONENT IN THE GRAPH, GENRCM OBTAINS THE ORDERING*/
13: /* BY CALLING THE SUBROUTINE RCM.*/
15: /* INPUT PARAMETERS -*/
16: /* NEQNS - NUMBER OF EQUATIONS*/
17: /* (XADJ, ADJNCY) - ARRAY PAIR CONTAINING THE ADJACENCY*/
18: /* STRUCTURE OF THE GRAPH OF THE MATRIX.*/
20: /* OUTPUT PARAMETER -*/
21: /* PERM - VECTOR THAT CONTAINS THE RCM ORDERING.*/
23: /* WORKING PARAMETERS -*/
24: /* MASK - IS USED TO MARK VARIABLES THAT HAVE BEEN*/
25: /* NUMBERED DURING THE ORDERING PROCESS. IT IS*/
26: /* INITIALIZED TO 1, AND SET TO ZERO AS EACH NODE*/
27: /* IS NUMBERED.*/
28: /* XLS - THE INDEX VECTOR FOR A LEVEL STRUCTURE. THE*/
29: /* LEVEL STRUCTURE IS STORED IN THE CURRENTLY*/
30: /* UNUSED SPACES IN THE PERMUTATION VECTOR PERM.*/
32: /* PROGRAM SUBROUTINES -*/
33: /* FN../../.., RCM.*/
34: /*****************************************************************/
37: PetscErrorCode SPARSEPACKgenrcm(PetscInt *neqns,PetscInt *xadj,PetscInt *adjncy,PetscInt *perm,PetscInt *mask,PetscInt *xls)
38: {
39: /* System generated locals */
40: PetscInt i__1;
42: /* Local variables */
43: PetscInt nlvl,root,i,ccsize;
44: EXTERN PetscErrorCode SPARSEPACKfnroot(PetscInt*,PetscInt *,PetscInt *,PetscInt *,PetscInt *,PetscInt *,PetscInt *),
45: SPARSEPACKrcm(PetscInt*,PetscInt *,PetscInt *,PetscInt *,PetscInt *,PetscInt *,PetscInt *);
46: PetscInt num;
49: /* Parameter adjustments */
50: --xls;
51: --mask;
52: --perm;
53: --adjncy;
54: --xadj;
56: i__1 = *neqns;
57: for (i = 1; i <= i__1; ++i) {
58: mask[i] = 1;
59: }
60: num = 1;
61: i__1 = *neqns;
62: for (i = 1; i <= i__1; ++i) {
63: /* FOR EACH MASKED CONNECTED COMPONENT ...*/
64: if (!mask[i]) {
65: goto L200;
66: }
67: root = i;
68: /* FIRST FIND A PSEUDO-PERIPHERAL NODE ../../...*/
69: /* NOTE THAT THE LEVEL STRUCTURE FOUND BY*/
70: /* FN../../.. IS STORED STARTING AT PERM(NUM).*/
71: /* THEN RCM IS CALLED TO ORDER THE COMPONENT*/
72: /* USING ../../.. AS THE STARTING NODE.*/
73: SPARSEPACKfnroot(&root,&xadj[1],&adjncy[1],&mask[1],&nlvl,&xls[1],&perm[num]);
74: SPARSEPACKrcm(&root,&xadj[1],&adjncy[1],&mask[1],&perm[num],&ccsize,&xls[1]);
75: num += ccsize;
76: if (num > *neqns) {
77: return(0);
78: }
79: L200:
80: ;
81: }
82: return(0);
83: }