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