1*8be712e4SBarry Smith /* qmdqt.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
6*8be712e4SBarry Smith /***************************************************************/
7*8be712e4SBarry Smith /******** QMDQT ..... QUOT MIN DEG QUOT TRANSFORM ********/
8*8be712e4SBarry Smith /***************************************************************/
9*8be712e4SBarry Smith
10*8be712e4SBarry Smith /* PURPOSE - THIS SUBROUTINE PERFORMS THE QUOTIENT GRAPH */
11*8be712e4SBarry Smith /* TRANSFORMATION AFTER A NODE HAS BEEN ELIMINATED.*/
12*8be712e4SBarry Smith
13*8be712e4SBarry Smith /* INPUT PARAMETERS -*/
14*8be712e4SBarry Smith /* ROOT - THE NODE JUST ELIMINATED. IT BECOMES THE*/
15*8be712e4SBarry Smith /* REPRESENTATIVE OF THE NEW SUPERNODE.*/
16*8be712e4SBarry Smith /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE.*/
17*8be712e4SBarry Smith /* (RCHSZE, RCHSET) - THE REACHABLE SET OF ROOT IN THE*/
18*8be712e4SBarry Smith /* OLD QUOTIENT GRAPH.*/
19*8be712e4SBarry Smith /* NBRHD - THE NEIGHBORHOOD SET WHICH WILL BE MERGED*/
20*8be712e4SBarry Smith /* WITH ROOT TO FORM THE NEW SUPERNODE.*/
21*8be712e4SBarry Smith /* MARKER - THE MARKER VECTOR.*/
22*8be712e4SBarry Smith
23*8be712e4SBarry Smith /* UPDATED PARAMETER -*/
24*8be712e4SBarry Smith /* ADJNCY - BECOMES THE ADJNCY OF THE QUOTIENT GRAPH.*/
25*8be712e4SBarry Smith /***************************************************************/
SPARSEPACKqmdqt(const PetscInt * root,const PetscInt * xadj,const PetscInt * inadjncy,PetscInt * marker,PetscInt * rchsze,PetscInt * rchset,PetscInt * nbrhd)26*8be712e4SBarry Smith PetscErrorCode SPARSEPACKqmdqt(const PetscInt *root, const PetscInt *xadj, const PetscInt *inadjncy, PetscInt *marker, PetscInt *rchsze, PetscInt *rchset, PetscInt *nbrhd)
27*8be712e4SBarry Smith {
28*8be712e4SBarry Smith PetscInt *adjncy = (PetscInt *)inadjncy; /* Used as temporary and reset within this function */
29*8be712e4SBarry Smith /* System generated locals */
30*8be712e4SBarry Smith PetscInt i__1, i__2;
31*8be712e4SBarry Smith
32*8be712e4SBarry Smith /* Local variables */
33*8be712e4SBarry Smith PetscInt inhd, irch, node, ilink, j, nabor, jstop, jstrt;
34*8be712e4SBarry Smith
35*8be712e4SBarry Smith PetscFunctionBegin;
36*8be712e4SBarry Smith /* Parameter adjustments */
37*8be712e4SBarry Smith --nbrhd;
38*8be712e4SBarry Smith --rchset;
39*8be712e4SBarry Smith --marker;
40*8be712e4SBarry Smith --adjncy;
41*8be712e4SBarry Smith --xadj;
42*8be712e4SBarry Smith
43*8be712e4SBarry Smith irch = 0;
44*8be712e4SBarry Smith inhd = 0;
45*8be712e4SBarry Smith node = *root;
46*8be712e4SBarry Smith L100:
47*8be712e4SBarry Smith jstrt = xadj[node];
48*8be712e4SBarry Smith jstop = xadj[node + 1] - 2;
49*8be712e4SBarry Smith if (jstop < jstrt) goto L300;
50*8be712e4SBarry Smith
51*8be712e4SBarry Smith /* PLACE REACH NODES INTO THE ADJACENT LIST OF NODE*/
52*8be712e4SBarry Smith i__1 = jstop;
53*8be712e4SBarry Smith for (j = jstrt; j <= i__1; ++j) {
54*8be712e4SBarry Smith ++irch;
55*8be712e4SBarry Smith adjncy[j] = rchset[irch];
56*8be712e4SBarry Smith if (irch >= *rchsze) goto L400;
57*8be712e4SBarry Smith }
58*8be712e4SBarry Smith /* LINK TO OTHER SPACE PROVIDED BY THE NBRHD SET.*/
59*8be712e4SBarry Smith L300:
60*8be712e4SBarry Smith ilink = adjncy[jstop + 1];
61*8be712e4SBarry Smith node = -ilink;
62*8be712e4SBarry Smith if (ilink < 0) goto L100;
63*8be712e4SBarry Smith ++inhd;
64*8be712e4SBarry Smith node = nbrhd[inhd];
65*8be712e4SBarry Smith adjncy[jstop + 1] = -node;
66*8be712e4SBarry Smith goto L100;
67*8be712e4SBarry Smith /* ALL REACHABLE NODES HAVE BEEN SAVED. END THE ADJ LIST.*/
68*8be712e4SBarry Smith /* ADD ROOT TO THE NBR LIST OF EACH NODE IN THE REACH SET.*/
69*8be712e4SBarry Smith L400:
70*8be712e4SBarry Smith adjncy[j + 1] = 0;
71*8be712e4SBarry Smith i__1 = *rchsze;
72*8be712e4SBarry Smith for (irch = 1; irch <= i__1; ++irch) {
73*8be712e4SBarry Smith node = rchset[irch];
74*8be712e4SBarry Smith if (marker[node] < 0) goto L600;
75*8be712e4SBarry Smith
76*8be712e4SBarry Smith jstrt = xadj[node];
77*8be712e4SBarry Smith jstop = xadj[node + 1] - 1;
78*8be712e4SBarry Smith i__2 = jstop;
79*8be712e4SBarry Smith for (j = jstrt; j <= i__2; ++j) {
80*8be712e4SBarry Smith nabor = adjncy[j];
81*8be712e4SBarry Smith if (marker[nabor] >= 0) goto L500;
82*8be712e4SBarry Smith adjncy[j] = *root;
83*8be712e4SBarry Smith goto L600;
84*8be712e4SBarry Smith L500:;
85*8be712e4SBarry Smith }
86*8be712e4SBarry Smith L600:;
87*8be712e4SBarry Smith }
88*8be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS);
89*8be712e4SBarry Smith }
90