Actual source code: qmdqt.c
petsc-3.9.4 2018-09-11
2: /* qmdqt.f -- translated by f2c (version 19931217).*/
4: #include <petscsys.h>
5: #include <petsc/private/matorderimpl.h>
7: /***************************************************************/
8: /******** QMDQT ..... QUOT MIN DEG QUOT TRANSFORM ********/
9: /***************************************************************/
11: /* PURPOSE - THIS SUBROUTINE PERFORMS THE QUOTIENT GRAPH */
12: /* TRANSFORMATION AFTER A NODE HAS BEEN ELIMINATED.*/
14: /* INPUT PARAMETERS -*/
15: /* ../../.. - THE NODE JUST ELIMINATED. IT BECOMES THE*/
16: /* REPRESENTATIVE OF THE NEW SUPERNODE.*/
17: /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE.*/
18: /* (RCHSZE, RCHSET) - THE REACHABLE SET OF ../../.. IN THE*/
19: /* OLD QUOTIENT GRAPH.*/
20: /* NBRHD - THE NEIGHBORHOOD SET WHICH WILL BE MERGED*/
21: /* WITH ../../.. TO FORM THE NEW SUPERNODE.*/
22: /* MARKER - THE MARKER VECTOR.*/
24: /* UPDATED PARAMETER -*/
25: /* ADJNCY - BECOMES THE ADJNCY OF THE QUOTIENT GRAPH.*/
26: /***************************************************************/
27: PetscErrorCode SPARSEPACKqmdqt(const PetscInt *root,const PetscInt *xadj,const PetscInt *inadjncy,PetscInt *marker, PetscInt *rchsze, PetscInt *rchset, PetscInt *nbrhd)
28: {
29: PetscInt *adjncy = (PetscInt*)inadjncy; /* Used as temporary and reset within this function */
30: /* System generated locals */
31: PetscInt i__1, i__2;
33: /* Local variables */
34: PetscInt inhd, irch, node, ilink, 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) goto L300;
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) goto L400;
58: }
59: /* LINK TO OTHER SPACE PROVIDED BY THE NBRHD SET.*/
60: L300:
61: ilink = adjncy[jstop + 1];
62: node = -ilink;
63: if (ilink < 0) goto L100;
64: ++inhd;
65: node = nbrhd[inhd];
66: adjncy[jstop + 1] = -node;
67: goto L100;
68: /* ALL REACHABLE NODES HAVE BEEN SAVED. END THE ADJ LIST.*/
69: /* ADD ../../.. TO THE NBR LIST OF EACH NODE IN THE REACH SET.*/
70: L400:
71: adjncy[j + 1] = 0;
72: i__1 = *rchsze;
73: for (irch = 1; irch <= i__1; ++irch) {
74: node = rchset[irch];
75: if (marker[node] < 0) goto L600;
77: jstrt = xadj[node];
78: jstop = xadj[node + 1] - 1;
79: i__2 = jstop;
80: for (j = jstrt; j <= i__2; ++j) {
81: nabor = adjncy[j];
82: if (marker[nabor] >= 0) goto L500;
83: adjncy[j] = *root;
84: goto L600;
85: L500:
86: ;
87: }
88: L600:
89: ;
90: }
91: return(0);
92: }