18be712e4SBarry Smith /* fn1wd.f -- translated by f2c (version 19931217).*/
28be712e4SBarry Smith
38be712e4SBarry Smith #include <petsc/private/matorderimpl.h>
48be712e4SBarry Smith
58be712e4SBarry Smith /*****************************************************************/
68be712e4SBarry Smith /******** FN1WD ..... FIND ONE-WAY DISSECTORS *********/
78be712e4SBarry Smith /*****************************************************************/
88be712e4SBarry Smith /* PURPOSE - THIS SUBROUTINE FINDS ONE-WAY DISSECTORS OF */
98be712e4SBarry Smith /* A CONNECTED COMPONENT SPECIFIED BY MASK AND ROOT. */
108be712e4SBarry Smith /* */
118be712e4SBarry Smith /* INPUT PARAMETERS - */
128be712e4SBarry Smith /* ROOT - A NODE THAT DEFINES (ALONG WITH MASK) THE */
138be712e4SBarry Smith /* COMPONENT TO BE PROCESSED. */
148be712e4SBarry Smith /* (XADJ, ADJNCY) - THE ADJACENCY STRUCTURE. */
158be712e4SBarry Smith /* */
168be712e4SBarry Smith /* OUTPUT PARAMETERS - */
178be712e4SBarry Smith /* NSEP - NUMBER OF NODES IN THE ONE-WAY DISSECTORS. */
188be712e4SBarry Smith /* SEP - VECTOR CONTAINING THE DISSECTOR NODES. */
198be712e4SBarry Smith /* */
208be712e4SBarry Smith /* UPDATED PARAMETER - */
218be712e4SBarry Smith /* MASK - NODES IN THE DISSECTOR HAVE THEIR MASK VALUES */
228be712e4SBarry Smith /* SET TO ZERO. */
238be712e4SBarry Smith /* */
248be712e4SBarry Smith /* WORKING PARAMETERS- */
258be712e4SBarry Smith /* (XLS, LS) - LEVEL STRUCTURE USED BY THE ROUTINE FNROOT. */
268be712e4SBarry Smith /* */
278be712e4SBarry Smith /* PROGRAM SUBROUTINE - */
288be712e4SBarry Smith /* FNROOT. */
298be712e4SBarry Smith /*****************************************************************/
SPARSEPACKfn1wd(PetscInt * root,const PetscInt * inxadj,const PetscInt * adjncy,PetscInt * mask,PetscInt * nsep,PetscInt * sep,PetscInt * nlvl,PetscInt * xls,PetscInt * ls)308be712e4SBarry Smith PetscErrorCode SPARSEPACKfn1wd(PetscInt *root, const PetscInt *inxadj, const PetscInt *adjncy, PetscInt *mask, PetscInt *nsep, PetscInt *sep, PetscInt *nlvl, PetscInt *xls, PetscInt *ls)
318be712e4SBarry Smith {
328be712e4SBarry Smith PetscInt *xadj = (PetscInt *)inxadj; /* Used as temporary and reset */
338be712e4SBarry Smith /* System generated locals */
348be712e4SBarry Smith PetscInt i__1, i__2;
358be712e4SBarry Smith
368be712e4SBarry Smith /* Local variables */
378be712e4SBarry Smith PetscInt node, i, j, k;
388be712e4SBarry Smith PetscReal width, fnlvl;
398be712e4SBarry Smith PetscInt kstop, kstrt, lp1beg, lp1end;
408be712e4SBarry Smith PetscReal deltp1;
418be712e4SBarry Smith PetscInt lvlbeg, lvlend;
428be712e4SBarry Smith PetscInt nbr, lvl;
438be712e4SBarry Smith
448be712e4SBarry Smith PetscFunctionBegin;
458be712e4SBarry Smith /* Parameter adjustments */
468be712e4SBarry Smith --ls;
478be712e4SBarry Smith --xls;
488be712e4SBarry Smith --sep;
498be712e4SBarry Smith --mask;
508be712e4SBarry Smith --adjncy;
518be712e4SBarry Smith --xadj;
528be712e4SBarry Smith
538be712e4SBarry Smith PetscCall(SPARSEPACKfnroot(root, &xadj[1], &adjncy[1], &mask[1], nlvl, &xls[1], &ls[1]));
548be712e4SBarry Smith fnlvl = (PetscReal)(*nlvl);
558be712e4SBarry Smith *nsep = xls[*nlvl + 1] - 1;
568be712e4SBarry Smith width = (PetscReal)(*nsep) / fnlvl;
578be712e4SBarry Smith deltp1 = PetscSqrtReal((width * 3. + 13.) / 2.) + 1.;
588be712e4SBarry Smith if (*nsep >= 50 && deltp1 <= fnlvl * .5f) goto L300;
598be712e4SBarry Smith
608be712e4SBarry Smith /* THE COMPONENT IS TOO SMALL, OR THE LEVEL STRUCTURE */
618be712e4SBarry Smith /* IS VERY LONG AND NARROW. RETURN THE WHOLE COMPONENT.*/
628be712e4SBarry Smith i__1 = *nsep;
638be712e4SBarry Smith for (i = 1; i <= i__1; ++i) {
648be712e4SBarry Smith node = ls[i];
658be712e4SBarry Smith sep[i] = node;
668be712e4SBarry Smith mask[node] = 0;
678be712e4SBarry Smith }
688be712e4SBarry Smith PetscFunctionReturn(PETSC_SUCCESS);
698be712e4SBarry Smith /* FIND THE PARALLEL DISSECTORS.*/
708be712e4SBarry Smith L300:
718be712e4SBarry Smith *nsep = 0;
728be712e4SBarry Smith i = 0;
738be712e4SBarry Smith L400:
748be712e4SBarry Smith ++i;
758be712e4SBarry Smith lvl = (PetscInt)((PetscReal)i * deltp1 + .5f);
768be712e4SBarry Smith if (lvl >= *nlvl) PetscFunctionReturn(PETSC_SUCCESS);
778be712e4SBarry Smith lvlbeg = xls[lvl];
788be712e4SBarry Smith lp1beg = xls[lvl + 1];
798be712e4SBarry Smith lvlend = lp1beg - 1;
808be712e4SBarry Smith lp1end = xls[lvl + 2] - 1;
818be712e4SBarry Smith i__1 = lp1end;
828be712e4SBarry Smith for (j = lp1beg; j <= i__1; ++j) {
838be712e4SBarry Smith node = ls[j];
848be712e4SBarry Smith xadj[node] = -xadj[node];
858be712e4SBarry Smith }
868be712e4SBarry Smith /* NODES IN LEVEL LVL ARE CHOSEN TO FORM DISSECTOR. */
878be712e4SBarry Smith /* INCLUDE ONLY THOSE WITH NEIGHBORS IN LVL+1 LEVEL. */
888be712e4SBarry Smith /* XADJ IS USED TEMPORARILY TO MARK NODES IN LVL+1. */
898be712e4SBarry Smith i__1 = lvlend;
908be712e4SBarry Smith for (j = lvlbeg; j <= i__1; ++j) {
918be712e4SBarry Smith node = ls[j];
928be712e4SBarry Smith kstrt = xadj[node];
938be712e4SBarry Smith i__2 = xadj[node + 1];
94*835f2295SStefano Zampini kstop = PetscAbsInt(i__2) - 1;
958be712e4SBarry Smith i__2 = kstop;
968be712e4SBarry Smith for (k = kstrt; k <= i__2; ++k) {
978be712e4SBarry Smith nbr = adjncy[k];
988be712e4SBarry Smith if (xadj[nbr] > 0) goto L600;
998be712e4SBarry Smith ++(*nsep);
1008be712e4SBarry Smith sep[*nsep] = node;
1018be712e4SBarry Smith mask[node] = 0;
1028be712e4SBarry Smith goto L700;
1038be712e4SBarry Smith L600:;
1048be712e4SBarry Smith }
1058be712e4SBarry Smith L700:;
1068be712e4SBarry Smith }
1078be712e4SBarry Smith i__1 = lp1end;
1088be712e4SBarry Smith for (j = lp1beg; j <= i__1; ++j) {
1098be712e4SBarry Smith node = ls[j];
1108be712e4SBarry Smith xadj[node] = -xadj[node];
1118be712e4SBarry Smith }
1128be712e4SBarry Smith goto L400;
1138be712e4SBarry Smith }
114