Actual source code: qmdqt.c
1: /* qmdqt.f -- translated by f2c (version 19931217).*/
3: #include petsc.h
5: /***************************************************************/
6: /******** QMDQT ..... QUOT MIN DEG QUOT TRANSFORM ********/
7: /***************************************************************/
9: /* PURPOSE - THIS SUBROUTINE PERFORMS THE QUOTIENT GRAPH */
10: /* TRANSFORMATION AFTER A NODE HAS BEEN ELIMINATED.*/
12: /* INPUT PARAMETERS -*/
13: /* ../../.. - THE NODE JUST ELIMINATED. IT BECOMES THE*/
14: /* REPRESENTATIVE OF THE NEW SUPERNODE.*/
15: /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE.*/
16: /* (RCHSZE, RCHSET) - THE REACHABLE SET OF ../../.. IN THE*/
17: /* OLD QUOTIENT GRAPH.*/
18: /* NBRHD - THE NEIGHBORHOOD SET WHICH WILL BE MERGED*/
19: /* WITH ../../.. TO FORM THE NEW SUPERNODE.*/
20: /* MARKER - THE MARKER VECTOR.*/
22: /* UPDATED PARAMETER -*/
23: /* ADJNCY - BECOMES THE ADJNCY OF THE QUOTIENT GRAPH.*/
24: /***************************************************************/
27: PetscErrorCode SPARSEPACKqmdqt(PetscInt *root, PetscInt *xadj, PetscInt *adjncy,
28: PetscInt *marker, PetscInt *rchsze, PetscInt *rchset, PetscInt *nbrhd)
29: {
30: /* System generated locals */
31: PetscInt i__1, i__2;
33: /* Local variables */
34: PetscInt inhd, irch, node, link, j, nabor, jstop, jstrt;
37: /* Parameter adjustments */
38: --nbrhd;
39: --rchset;
40: --marker;
41: --adjncy;
42: --xadj;
44: irch = 0;
45: inhd = 0;
46: node = *root;
47: L100:
48: jstrt = xadj[node];
49: jstop = xadj[node + 1] - 2;
50: if (jstop < jstrt) {
51: goto L300;
52: }
53: /* PLACE REACH NODES INTO THE ADJACENT LIST OF NODE*/
54: i__1 = jstop;
55: for (j = jstrt; j <= i__1; ++j) {
56: ++irch;
57: adjncy[j] = rchset[irch];
58: if (irch >= *rchsze) {
59: goto L400;
60: }
61: }
62: /* LINK TO OTHER SPACE PROVIDED BY THE NBRHD SET.*/
63: L300:
64: link = adjncy[jstop + 1];
65: node = -link;
66: if (link < 0) {
67: goto L100;
68: }
69: ++inhd;
70: node = nbrhd[inhd];
71: adjncy[jstop + 1] = -node;
72: goto L100;
73: /* ALL REACHABLE NODES HAVE BEEN SAVED. END THE ADJ LIST.*/
74: /* ADD ../../.. TO THE NBR LIST OF EACH NODE IN THE REACH SET.*/
75: L400:
76: adjncy[j + 1] = 0;
77: i__1 = *rchsze;
78: for (irch = 1; irch <= i__1; ++irch) {
79: node = rchset[irch];
80: if (marker[node] < 0) {
81: goto L600;
82: }
83: jstrt = xadj[node];
84: jstop = xadj[node + 1] - 1;
85: i__2 = jstop;
86: for (j = jstrt; j <= i__2; ++j) {
87: nabor = adjncy[j];
88: if (marker[nabor] >= 0) {
89: goto L500;
90: }
91: adjncy[j] = *root;
92: goto L600;
93: L500:
94: ;
95: }
96: L600:
97: ;
98: }
99: return(0);
100: }