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