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