18be712e4SBarry Smith /* degree.f -- translated by f2c (version 19931217).*/
28be712e4SBarry Smith
38be712e4SBarry Smith #include <petsc/private/matorderimpl.h>
48be712e4SBarry Smith
58be712e4SBarry Smith /*****************************************************************/
68be712e4SBarry Smith /********* DEGREE ..... DEGREE IN MASKED COMPONENT *********/
78be712e4SBarry Smith /*****************************************************************/
88be712e4SBarry Smith
98be712e4SBarry Smith /* PURPOSE - THIS ROUTINE COMPUTES THE DEGREES OF THE NODES*/
108be712e4SBarry Smith /* IN THE CONNECTED COMPONENT SPECIFIED BY MASK AND ROOT*/
118be712e4SBarry Smith /* NODES FOR WHICH MASK IS ZERO ARE IGNORED.*/
128be712e4SBarry Smith
138be712e4SBarry Smith /* INPUT PARAMETER -*/
148be712e4SBarry Smith /* ROOT - IS THE INPUT NODE THAT DEFINES THE COMPONENT.*/
158be712e4SBarry Smith /* (XADJ, ADJNCY) - ADJACENCY STRUCTURE PAIR.*/
168be712e4SBarry Smith /* MASK - SPECIFIES A SECTION SUBGRAPH.*/
178be712e4SBarry Smith
188be712e4SBarry Smith /* OUTPUT PARAMETERS -*/
198be712e4SBarry Smith /* DEG - ARRAY CONTAINING THE DEGREES OF THE NODES IN*/
208be712e4SBarry Smith /* THE COMPONENT.*/
218be712e4SBarry Smith /* CCSIZE-SIZE OF THE COMPONENT SPECIFIED BY MASK AND ROOT*/
228be712e4SBarry Smith /* WORKING PARAMETER -*/
238be712e4SBarry Smith /* LS - A TEMPORARY VECTOR USED TO STORE THE NODES OF THE*/
248be712e4SBarry Smith /* COMPONENT LEVEL BY LEVEL.*/
258be712e4SBarry Smith /*****************************************************************/
SPARSEPACKdegree(const PetscInt * root,const PetscInt * inxadj,const PetscInt * adjncy,PetscInt * mask,PetscInt * deg,PetscInt * ccsize,PetscInt * ls)268be712e4SBarry Smith PetscErrorCode SPARSEPACKdegree(const PetscInt *root, const PetscInt *inxadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *deg, PetscInt *ccsize, PetscInt *ls)
278be712e4SBarry Smith {
288be712e4SBarry Smith PetscInt *xadj = (PetscInt *)inxadj; /* Used as temporary and reset within this function */
298be712e4SBarry Smith /* System generated locals */
308be712e4SBarry Smith PetscInt i__1, i__2;
318be712e4SBarry Smith
328be712e4SBarry Smith /* Local variables */
338be712e4SBarry Smith PetscInt ideg, node, i, j, jstop, jstrt, lbegin, lvlend, lvsize, nbr;
348be712e4SBarry Smith /* INITIALIZATION ...*/
358be712e4SBarry Smith /* THE ARRAY XADJ IS USED AS A TEMPORARY MARKER TO*/
368be712e4SBarry Smith /* INDICATE WHICH NODES HAVE BEEN CONSIDERED SO FAR.*/
378be712e4SBarry Smith
388be712e4SBarry Smith PetscFunctionBegin;
398be712e4SBarry Smith /* Parameter adjustments */
408be712e4SBarry Smith --ls;
418be712e4SBarry Smith --deg;
428be712e4SBarry Smith --mask;
438be712e4SBarry Smith --adjncy;
448be712e4SBarry Smith --xadj;
458be712e4SBarry Smith
468be712e4SBarry Smith ls[1] = *root;
478be712e4SBarry Smith xadj[*root] = -xadj[*root];
488be712e4SBarry Smith lvlend = 0;
498be712e4SBarry Smith *ccsize = 1;
508be712e4SBarry Smith /* LBEGIN IS THE POINTER TO THE BEGINNING OF THE CURRENT*/
518be712e4SBarry Smith /* LEVEL, AND LVLEND POINTS TO THE END OF THIS LEVEL.*/
528be712e4SBarry Smith L100:
538be712e4SBarry Smith lbegin = lvlend + 1;
548be712e4SBarry Smith lvlend = *ccsize;
558be712e4SBarry Smith /* FIND THE DEGREES OF NODES IN THE CURRENT LEVEL,*/
568be712e4SBarry Smith /* AND AT THE SAME TIME, GENERATE THE NEXT LEVEL.*/
578be712e4SBarry Smith i__1 = lvlend;
588be712e4SBarry Smith for (i = lbegin; i <= i__1; ++i) {
598be712e4SBarry Smith node = ls[i];
608be712e4SBarry Smith jstrt = -xadj[node];
618be712e4SBarry Smith i__2 = xadj[node + 1];
62*835f2295SStefano Zampini jstop = PetscAbsInt(i__2) - 1;
638be712e4SBarry Smith ideg = 0;
648be712e4SBarry Smith if (jstop < jstrt) goto L300;
658be712e4SBarry Smith i__2 = jstop;
668be712e4SBarry Smith for (j = jstrt; j <= i__2; ++j) {
678be712e4SBarry Smith nbr = adjncy[j];
688be712e4SBarry Smith if (!mask[nbr]) goto L200;
698be712e4SBarry Smith ++ideg;
708be712e4SBarry Smith if (xadj[nbr] < 0) goto L200;
718be712e4SBarry Smith xadj[nbr] = -xadj[nbr];
728be712e4SBarry Smith ++(*ccsize);
738be712e4SBarry Smith ls[*ccsize] = nbr;
748be712e4SBarry Smith L200:;
758be712e4SBarry Smith }
768be712e4SBarry Smith L300:
778be712e4SBarry Smith deg[node] = ideg;
788be712e4SBarry Smith }
798be712e4SBarry Smith /* COMPUTE THE CURRENT LEVEL WIDTH. */
808be712e4SBarry Smith /* IF IT IS NONZERO, GENERATE ANOTHER LEVEL.*/
818be712e4SBarry Smith lvsize = *ccsize - lvlend;
828be712e4SBarry Smith if (lvsize > 0) goto L100;
838be712e4SBarry Smith /* RESET XADJ TO ITS CORRECT SIGN AND RETURN. */
848be712e4SBarry Smith i__1 = *ccsize;
858be712e4SBarry Smith for (i = 1; i <= i__1; ++i) {
868be712e4SBarry Smith node = ls[i];
878be712e4SBarry Smith xadj[node] = -xadj[node];
888be712e4SBarry Smith }
898be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS);
908be712e4SBarry Smith }
91