C++ ComputeCut函数代码示例

/*************************************************************************
* This function performs k-way refinement
**************************************************************************/
void MCGreedy_KWayEdgeBalanceHorizontal(CtrlType *ctrl, GraphType *graph, int nparts, 
       float *ubvec, int npasses)
{
  int i, ii, /*iii,*/ j, /*jj,*/ k, /*l,*/ pass, nvtxs, ncon, nbnd, myndegrees, oldgain, gain, nmoves; 
  int from, me, to, oldcut;
  idxtype *xadj, *adjncy, *adjwgt;
  idxtype *where, *perm, *bndptr, *bndind, *moved;
  EDegreeType *myedegrees;
  RInfoType *myrinfo;
  PQueueType queue;
  float *npwgts, *nvwgt, *minwgt, *maxwgt, tvec[MAXNCON];

  nvtxs = graph->nvtxs;
  ncon = graph->ncon;
  xadj = graph->xadj;
  adjncy = graph->adjncy;
  adjwgt = graph->adjwgt;

  bndind = graph->bndind;
  bndptr = graph->bndptr;

  where = graph->where;
  npwgts = graph->npwgts;
  
  /* Setup the weight intervals of the various subdomains */
  minwgt =  fwspacemalloc(ctrl, ncon*nparts);
  maxwgt = fwspacemalloc(ctrl, ncon*nparts);

  for (i=0; i<nparts; i++) {
    for (j=0; j<ncon; j++) {
      maxwgt[i*ncon+j] = ubvec[j]/nparts;
      minwgt[i*ncon+j] = 1.0/(ubvec[j]*nparts);
    }
  }

  perm = idxwspacemalloc(ctrl, nvtxs);
  moved = idxwspacemalloc(ctrl, nvtxs);

  PQueueInit(ctrl, &queue, nvtxs, graph->adjwgtsum[idxamax(nvtxs, graph->adjwgtsum)]);

  if (ctrl->dbglvl&DBG_REFINE) {
    printf("Partitions: [%5.4f %5.4f], Nv-Nb[%6d %6d]. Cut: %6d, LB: ",
            npwgts[samin(ncon*nparts, npwgts)], npwgts[samax(ncon*nparts, npwgts)], 
            graph->nvtxs, graph->nbnd, graph->mincut);
    ComputeHKWayLoadImbalance(ncon, nparts, npwgts, tvec);
    for (i=0; i<ncon; i++)
      printf("%.3f ", tvec[i]);
    printf("[B]\n");
  }


  for (pass=0; pass<npasses; pass++) {
    ASSERT(ComputeCut(graph, where) == graph->mincut);

    /* Check to see if things are out of balance, given the tolerance */
    if (MocIsHBalanced(ncon, nparts, npwgts, ubvec))
      break;

    PQueueReset(&queue);
    idxset(nvtxs, -1, moved);

    oldcut = graph->mincut;
    nbnd = graph->nbnd;

    RandomPermute(nbnd, perm, 1);
    for (ii=0; ii<nbnd; ii++) {
      i = bndind[perm[ii]];
      PQueueInsert(&queue, i, graph->rinfo[i].ed - graph->rinfo[i].id);
      moved[i] = 2;
    }

    nmoves = 0;
    for (;;) {
      if ((i = PQueueGetMax(&queue)) == -1) 
        break;
      moved[i] = 1;

      myrinfo = graph->rinfo+i;
      from = where[i];
      nvwgt = graph->nvwgt+i*ncon;

      if (AreAllHVwgtsBelow(ncon, 1.0, npwgts+from*ncon, -1.0, nvwgt, minwgt+from*ncon))
        continue;   /* This cannot be moved! */

      myedegrees = myrinfo->edegrees;
      myndegrees = myrinfo->ndegrees;

      for (k=0; k<myndegrees; k++) {
        to = myedegrees[k].pid;
        if (IsHBalanceBetterFT(ncon, nparts, npwgts+from*ncon, npwgts+to*ncon, nvwgt, ubvec))
          break;
      }
      if (k == myndegrees) 
        continue;  /* break out if you did not find a candidate */

      for (j=k+1; j<myndegrees; j++) {
        to = myedegrees[j].pid;
//.........这里部分代码省略.........

void FM_2WayCutRefine(ctrl_t *ctrl, graph_t *graph, real_t *ntpwgts, idx_t niter)
{
    idx_t i, ii, j, k, kwgt, nvtxs, nbnd, nswaps, from, to, pass, me, limit, tmp;
    idx_t *xadj, *vwgt, *adjncy, *adjwgt, *where, *id, *ed, *bndptr, *bndind, *pwgts;
    idx_t *moved, *swaps, *perm;
    rpq_t *queues[2];
    idx_t higain, mincut, mindiff, origdiff, initcut, newcut, mincutorder, avgvwgt;
    idx_t tpwgts[2];

    WCOREPUSH;

    nvtxs  = graph->nvtxs;
    xadj   = graph->xadj;
    vwgt   = graph->vwgt;
    adjncy = graph->adjncy;
    adjwgt = graph->adjwgt;
    where  = graph->where;
    id     = graph->id;
    ed     = graph->ed;
    pwgts  = graph->pwgts;
    bndptr = graph->bndptr;
    bndind = graph->bndind;

    moved = iwspacemalloc(ctrl, nvtxs);
    swaps = iwspacemalloc(ctrl, nvtxs);
    perm  = iwspacemalloc(ctrl, nvtxs);

    tpwgts[0] = graph->tvwgt[0]*ntpwgts[0];
    tpwgts[1] = graph->tvwgt[0]-tpwgts[0];

    limit   = gk_min(gk_max(0.01*nvtxs, 15), 100);
    avgvwgt = gk_min((pwgts[0]+pwgts[1])/20, 2*(pwgts[0]+pwgts[1])/nvtxs);

    queues[0] = rpqCreate(nvtxs);
    queues[1] = rpqCreate(nvtxs);

    IFSET(ctrl->dbglvl, METIS_DBG_REFINE,
          Print2WayRefineStats(ctrl, graph, ntpwgts, 0, -2));

    origdiff = iabs(tpwgts[0]-pwgts[0]);
    iset(nvtxs, -1, moved);
    for (pass=0; pass<niter; pass++) { /* Do a number of passes */
        rpqReset(queues[0]);
        rpqReset(queues[1]);

        mincutorder = -1;
        newcut = mincut = initcut = graph->mincut;
        mindiff = iabs(tpwgts[0]-pwgts[0]);

        ASSERT(ComputeCut(graph, where) == graph->mincut);
        ASSERT(CheckBnd(graph));

        /* Insert boundary nodes in the priority queues */
        nbnd = graph->nbnd;
        irandArrayPermute(nbnd, perm, nbnd, 1);
        for (ii=0; ii<nbnd; ii++) {
            i = perm[ii];
            ASSERT(ed[bndind[i]] > 0 || id[bndind[i]] == 0);
            ASSERT(bndptr[bndind[i]] != -1);
            rpqInsert(queues[where[bndind[i]]], bndind[i], ed[bndind[i]]-id[bndind[i]]);
        }

        for (nswaps=0; nswaps<nvtxs; nswaps++) {
            from = (tpwgts[0]-pwgts[0] < tpwgts[1]-pwgts[1] ? 0 : 1);
            to = (from+1)%2;

            if ((higain = rpqGetTop(queues[from])) == -1)
                break;
            ASSERT(bndptr[higain] != -1);

            newcut -= (ed[higain]-id[higain]);
            INC_DEC(pwgts[to], pwgts[from], vwgt[higain]);

            if ((newcut < mincut && iabs(tpwgts[0]-pwgts[0]) <= origdiff+avgvwgt) ||
                    (newcut == mincut && iabs(tpwgts[0]-pwgts[0]) < mindiff)) {
                mincut  = newcut;
                mindiff = iabs(tpwgts[0]-pwgts[0]);
                mincutorder = nswaps;
            }
            else if (nswaps-mincutorder > limit) { /* We hit the limit, undo last move */
                newcut += (ed[higain]-id[higain]);
                INC_DEC(pwgts[from], pwgts[to], vwgt[higain]);
                break;
            }

            where[higain] = to;
            moved[higain] = nswaps;
            swaps[nswaps] = higain;

            IFSET(ctrl->dbglvl, METIS_DBG_MOVEINFO,
                  printf("Moved %6"PRIDX" from %"PRIDX". [%3"PRIDX" %3"PRIDX"] %5"PRIDX" [%4"PRIDX" %4"PRIDX"]\n", higain, from, ed[higain]-id[higain], vwgt[higain], newcut, pwgts[0], pwgts[1]));

            /**************************************************************
            * Update the id[i]/ed[i] values of the affected nodes
            ***************************************************************/
            SWAP(id[higain], ed[higain], tmp);
            if (ed[higain] == 0 && xadj[higain] < xadj[higain+1])
                BNDDelete(nbnd, bndind,  bndptr, higain);

            for (j=xadj[higain]; j<xadj[higain+1]; j++) {
//.........这里部分代码省略.........

/*************************************************************************
* This function performs k-way refinement
**************************************************************************/
void MCRandom_KWayEdgeRefineHorizontal(CtrlType *ctrl, GraphType *graph, int nparts, 
       float *orgubvec, int npasses)
{
  int i, ii, iii, j, /*jj,*/ k, /*l,*/ pass, nvtxs, ncon, nmoves, nbnd, myndegrees, same; 
  int from, me, to, oldcut, gain;
  idxtype *xadj, *adjncy, *adjwgt;
  idxtype *where, *perm, *bndptr, *bndind;
  EDegreeType *myedegrees;
  RInfoType *myrinfo;
  float *npwgts, *nvwgt, *minwgt, *maxwgt, maxlb, minlb, ubvec[MAXNCON], tvec[MAXNCON];

  nvtxs = graph->nvtxs;
  ncon = graph->ncon;
  xadj = graph->xadj;
  adjncy = graph->adjncy;
  adjwgt = graph->adjwgt;

  bndptr = graph->bndptr;
  bndind = graph->bndind;

  where = graph->where;
  npwgts = graph->npwgts;
  
  /* Setup the weight intervals of the various subdomains */
  minwgt =  fwspacemalloc(ctrl, nparts*ncon);
  maxwgt = fwspacemalloc(ctrl, nparts*ncon);

  /* See if the orgubvec consists of identical constraints */
  maxlb = minlb = orgubvec[0];
  for (i=1; i<ncon; i++) {
    minlb = (orgubvec[i] < minlb ? orgubvec[i] : minlb);
    maxlb = (orgubvec[i] > maxlb ? orgubvec[i] : maxlb);
  }
  same = (fabs(maxlb-minlb) < .01 ? 1 : 0);


  /* Let's not get very optimistic. Let Balancing do the work */
  ComputeHKWayLoadImbalance(ncon, nparts, npwgts, ubvec);
  for (i=0; i<ncon; i++)
    ubvec[i] = amax(ubvec[i], orgubvec[i]);

  if (!same) {
    for (i=0; i<nparts; i++) {
      for (j=0; j<ncon; j++) {
        maxwgt[i*ncon+j] = ubvec[j]/nparts;
        minwgt[i*ncon+j] = 1.0/(ubvec[j]*nparts);
      }
    }
  }
  else {
    maxlb = ubvec[0];
    for (i=1; i<ncon; i++) 
      maxlb = (ubvec[i] > maxlb ? ubvec[i] : maxlb);

    for (i=0; i<nparts; i++) {
      for (j=0; j<ncon; j++) {
        maxwgt[i*ncon+j] = maxlb/nparts;
        minwgt[i*ncon+j] = 1.0/(maxlb*nparts);
      }
    }
  }


  perm = idxwspacemalloc(ctrl, nvtxs);

  if (ctrl->dbglvl&DBG_REFINE) {
    printf("Partitions: [%5.4f %5.4f], Nv-Nb[%6d %6d]. Cut: %6d, LB: ",
            npwgts[samin(ncon*nparts, npwgts)], npwgts[samax(ncon*nparts, npwgts)], 
            graph->nvtxs, graph->nbnd, graph->mincut);
    ComputeHKWayLoadImbalance(ncon, nparts, npwgts, tvec);
    for (i=0; i<ncon; i++)
      printf("%.3f ", tvec[i]);
    printf("\n");
  }

  for (pass=0; pass<npasses; pass++) {
    ASSERT(ComputeCut(graph, where) == graph->mincut);

    oldcut = graph->mincut;
    nbnd = graph->nbnd;

    RandomPermute(nbnd, perm, 1);
    for (nmoves=iii=0; iii<graph->nbnd; iii++) {
      ii = perm[iii];
      if (ii >= nbnd)
        continue;
      i = bndind[ii];

      myrinfo = graph->rinfo+i;

      if (myrinfo->ed >= myrinfo->id) { /* Total ED is too high */
        from = where[i];
        nvwgt = graph->nvwgt+i*ncon;

        if (myrinfo->id > 0 && AreAllHVwgtsBelow(ncon, 1.0, npwgts+from*ncon, -1.0, nvwgt, minwgt+from*ncon)) 
          continue;   /* This cannot be moved! */

//.........这里部分代码省略.........

/*************************************************************************
* This function performs an edge-based FM refinement
**************************************************************************/
void MocGeneral2WayBalance2(CtrlType *ctrl, GraphType *graph, float *tpwgts, float *ubvec)
{
  int i, ii, j, k, l, kwgt, nvtxs, ncon, nbnd, nswaps, from, to, limit, tmp, cnum;
  idxtype *xadj, *adjncy, *adjwgt, *where, *id, *ed, *bndptr, *bndind;
  idxtype *moved, *swaps, *perm, *qnum;
  float *nvwgt, *npwgts, origbal[MAXNCON], minbal[MAXNCON], newbal[MAXNCON];
  PQueueType parts[MAXNCON][2];
  int higain, oldgain, mincut, newcut, mincutorder;
  float *maxwgt, *minwgt, tvec[MAXNCON];


  nvtxs = graph->nvtxs;
  ncon = graph->ncon;
  xadj = graph->xadj;
  nvwgt = graph->nvwgt;
  adjncy = graph->adjncy;
  adjwgt = graph->adjwgt;
  where = graph->where;
  id = graph->id;
  ed = graph->ed;
  npwgts = graph->npwgts;
  bndptr = graph->bndptr;
  bndind = graph->bndind;

  moved = idxwspacemalloc(ctrl, nvtxs);
  swaps = idxwspacemalloc(ctrl, nvtxs);
  perm = idxwspacemalloc(ctrl, nvtxs);
  qnum = idxwspacemalloc(ctrl, nvtxs);

  limit = amin(amax(0.01*nvtxs, 15), 100);

  /* Setup the weight intervals of the two subdomains */
  minwgt = fwspacemalloc(ctrl, 2*ncon);
  maxwgt = fwspacemalloc(ctrl, 2*ncon);

  for (i=0; i<2; i++) {
    for (j=0; j<ncon; j++) {
      maxwgt[i*ncon+j] = tpwgts[i]*ubvec[j];
      minwgt[i*ncon+j] = tpwgts[i]*(1.0/ubvec[j]);
    }
  }


  /* Initialize the queues */
  for (i=0; i<ncon; i++) {
    PQueueInit(ctrl, &parts[i][0], nvtxs, PLUS_GAINSPAN+1);
    PQueueInit(ctrl, &parts[i][1], nvtxs, PLUS_GAINSPAN+1);
  }
  for (i=0; i<nvtxs; i++)
    qnum[i] = samax(ncon, nvwgt+i*ncon);

  Compute2WayHLoadImbalanceVec(ncon, npwgts, tpwgts, origbal);
  for (i=0; i<ncon; i++) 
    minbal[i] = origbal[i];

  newcut = mincut = graph->mincut;
  mincutorder = -1;

  if (ctrl->dbglvl&DBG_REFINE) {
    printf("Parts: [");
    for (l=0; l<ncon; l++)
      printf("(%.3f, %.3f) ", npwgts[l], npwgts[ncon+l]);
    printf("] T[%.3f %.3f], Nv-Nb[%5d, %5d]. ICut: %6d, LB: ", tpwgts[0], tpwgts[1], 
            graph->nvtxs, graph->nbnd, graph->mincut);
    for (i=0; i<ncon; i++)
      printf("%.3f ", origbal[i]);
    printf("[B]\n");
  }

  idxset(nvtxs, -1, moved);

  ASSERT(ComputeCut(graph, where) == graph->mincut);
  ASSERT(CheckBnd(graph));

  /* Insert all nodes in the priority queues */
  nbnd = graph->nbnd;
  RandomPermute(nvtxs, perm, 1);
  for (ii=0; ii<nvtxs; ii++) {
    i = perm[ii];
    PQueueInsert(&parts[qnum[i]][where[i]], i, ed[i]-id[i]);
  }


  for (nswaps=0; nswaps<nvtxs; nswaps++) {
    if (AreAllBelow(ncon, minbal, ubvec))
      break;

    SelectQueue3(ncon, npwgts, tpwgts, &from, &cnum, parts, maxwgt);
    to = (from+1)%2;

    if (from == -1 || (higain = PQueueGetMax(&parts[cnum][from])) == -1)
      break;

    saxpy(ncon, 1.0, nvwgt+higain*ncon, 1, npwgts+to*ncon, 1);
    saxpy(ncon, -1.0, nvwgt+higain*ncon, 1, npwgts+from*ncon, 1);
    newcut -= (ed[higain]-id[higain]);
    Compute2WayHLoadImbalanceVec(ncon, npwgts, tpwgts, newbal);
//.........这里部分代码省略.........

void FM_Mc2WayCutRefine(ctrl_t *ctrl, graph_t *graph, real_t *ntpwgts, idx_t niter)
{
    idx_t i, ii, j, k, l, kwgt, nvtxs, ncon, nbnd, nswaps, from, to, pass,
          me, limit, tmp, cnum;
    idx_t *xadj, *adjncy, *vwgt, *adjwgt, *pwgts, *where, *id, *ed,
          *bndptr, *bndind;
    idx_t *moved, *swaps, *perm, *qnum;
    idx_t higain, mincut, initcut, newcut, mincutorder;
    real_t *invtvwgt, *ubfactors, *minbalv, *newbalv;
    real_t origbal, minbal, newbal, rgain, ffactor;
    rpq_t **queues;

    WCOREPUSH;

    nvtxs    = graph->nvtxs;
    ncon     = graph->ncon;
    xadj     = graph->xadj;
    vwgt     = graph->vwgt;
    adjncy   = graph->adjncy;
    adjwgt   = graph->adjwgt;
    invtvwgt = graph->invtvwgt;
    where    = graph->where;
    id       = graph->id;
    ed       = graph->ed;
    pwgts    = graph->pwgts;
    bndptr   = graph->bndptr;
    bndind   = graph->bndind;

    moved     = iwspacemalloc(ctrl, nvtxs);
    swaps     = iwspacemalloc(ctrl, nvtxs);
    perm      = iwspacemalloc(ctrl, nvtxs);
    qnum      = iwspacemalloc(ctrl, nvtxs);
    ubfactors = rwspacemalloc(ctrl, ncon);
    newbalv   = rwspacemalloc(ctrl, ncon);
    minbalv   = rwspacemalloc(ctrl, ncon);

    limit = gk_min(gk_max(0.01*nvtxs, 25), 150);


    /* Determine a fudge factor to allow the refinement routines to get out
       of tight balancing constraints. */
    ffactor = .5/gk_max(20, nvtxs);

    /* Initialize the queues */
    queues = (rpq_t **)wspacemalloc(ctrl, 2*ncon*sizeof(rpq_t *));
    for (i=0; i<2*ncon; i++)
        queues[i] = rpqCreate(nvtxs);
    for (i=0; i<nvtxs; i++)
        qnum[i] = iargmax_nrm(ncon, vwgt+i*ncon, invtvwgt);

    /* Determine the unbalance tolerance for each constraint. The tolerance is
       equal to the maximum of the original load imbalance and the user-supplied
       allowed tolerance. The rationale behind this approach is to allow the
       refinement routine to improve the cut, without having to worry about fixing
       load imbalance problems. The load imbalance is addressed by the balancing
       routines. */
    origbal = ComputeLoadImbalanceDiffVec(graph, 2, ctrl->pijbm, ctrl->ubfactors, ubfactors);
    for (i=0; i<ncon; i++)
        ubfactors[i] = (ubfactors[i] > 0 ? ctrl->ubfactors[i]+ubfactors[i] : ctrl->ubfactors[i]);


    IFSET(ctrl->dbglvl, METIS_DBG_REFINE,
          Print2WayRefineStats(ctrl, graph, ntpwgts, origbal, -2));

    iset(nvtxs, -1, moved);
    for (pass=0; pass<niter; pass++) { /* Do a number of passes */
        for (i=0; i<2*ncon; i++)
            rpqReset(queues[i]);

        mincutorder = -1;
        newcut = mincut = initcut = graph->mincut;

        minbal = ComputeLoadImbalanceDiffVec(graph, 2, ctrl->pijbm, ubfactors, minbalv);

        ASSERT(ComputeCut(graph, where) == graph->mincut);
        ASSERT(CheckBnd(graph));

        /* Insert boundary nodes in the priority queues */
        nbnd = graph->nbnd;
        irandArrayPermute(nbnd, perm, nbnd/5, 1);
        for (ii=0; ii<nbnd; ii++) {
            i = bndind[perm[ii]];
            ASSERT(ed[i] > 0 || id[i] == 0);
            ASSERT(bndptr[i] != -1);
            //rgain = 1.0*(ed[i]-id[i])/sqrt(vwgt[i*ncon+qnum[i]]+1);
            //rgain = (ed[i]-id[i] > 0 ? 1.0*(ed[i]-id[i])/sqrt(vwgt[i*ncon+qnum[i]]+1) : ed[i]-id[i]);
            rgain = ed[i]-id[i];
            rpqInsert(queues[2*qnum[i]+where[i]], i, rgain);
        }

        for (nswaps=0; nswaps<nvtxs; nswaps++) {
            SelectQueue(graph, ctrl->pijbm, ubfactors, queues, &from, &cnum);

            to = (from+1)%2;

            if (from == -1 || (higain = rpqGetTop(queues[2*cnum+from])) == -1)
                break;
            ASSERT(bndptr[higain] != -1);

            newcut -= (ed[higain]-id[higain]);
//.........这里部分代码省略.........

/*************************************************************************
* This function is the entry point for KMETIS
**************************************************************************/
void *METIS_PartGraphForContact(idxtype *nvtxs, idxtype *xadj, idxtype *adjncy, 
                double *xyzcoords, idxtype *sflag, idxtype *numflag, idxtype *nparts, 
                idxtype *options, idxtype *edgecut, idxtype *part) 
{
  idxtype i, j, ii, dim, ncon, wgtflag, mcnumflag, nnodes, nlnodes, nclean, naclean, ndirty, maxdepth, rwgtflag, rnumflag;
  idxtype *mcvwgt, *dtpart, *marker, *leafpart;
  idxtype *adjwgt;
  float rubvec[2], lbvec[2];
  GraphType graph, *cgraph;
  ContactInfoType *cinfo;
  DKeyValueType *xyzcand[3];

  if (*numflag == 1)
    Change2CNumbering(*nvtxs, xadj, adjncy);

  /*---------------------------------------------------------------------
   * Allocate memory for the contact info type
   *---------------------------------------------------------------------*/
  cinfo = (ContactInfoType *)gk_malloc(sizeof(ContactInfoType), "METIS_PartGraphForContact: cinfo");
  cinfo->leafptr  = idxsmalloc(*nvtxs+1, 0, "METIS_PartGraphForContact: leafptr");
  cinfo->leafind  = idxsmalloc(*nvtxs, 0, "METIS_PartGraphForContact: leafind");
  cinfo->leafwgt  = idxsmalloc(*nvtxs, 0, "METIS_PartGraphForContact: leafwgt");
  cinfo->part     = idxsmalloc(*nvtxs, 0, "METIS_PartGraphForContact: part");
  leafpart = cinfo->leafpart = idxmalloc(*nvtxs, "METIS_PartGraphForContact: leafpart");
  cinfo->dtree    = (DTreeNodeType *)gk_malloc(sizeof(DTreeNodeType)*(*nvtxs), "METIS_PartGraphForContact: cinfo->dtree");
  cinfo->nvtxs    = *nvtxs;

  /*---------------------------------------------------------------------
   * Compute the initial k-way partitioning 
   *---------------------------------------------------------------------*/
  mcvwgt = idxsmalloc(2*(*nvtxs), 0, "METIS_PartGraphForContact: mcvwgt");
  for (i=0; i<*nvtxs; i++) {
    mcvwgt[2*i+0] = 1;
    mcvwgt[2*i+1] = (sflag[i] == 0 ? 0 : 1);
  }

  adjwgt = idxmalloc(xadj[*nvtxs], "METIS_PartGraphForContact: adjwgt");
  for (i=0; i<*nvtxs; i++) {
    for (j=xadj[i]; j<xadj[i+1]; j++) 
      adjwgt[j] = (sflag[i] && sflag[adjncy[j]] ? 5 : 1);
  }

  rubvec[0] = 1.03;
  rubvec[1] = 1.05;
  ncon = 2;
  mcnumflag = 0;
  wgtflag   = 1;

  METIS_mCPartGraphKway(nvtxs, &ncon, xadj, adjncy, mcvwgt, adjwgt, &wgtflag, &mcnumflag,
                        nparts, rubvec, options, edgecut, part);

  /* The following is just for stat reporting purposes */
  SetUpGraph(&graph, OP_KMETIS, *nvtxs, 2, xadj, adjncy, mcvwgt, NULL, 0);
  graph.vwgt = mcvwgt;
  ComputePartitionBalance(&graph, *nparts, part, lbvec);
  mprintf("  %D-way Edge-Cut: %7D, Balance: %5.2f %5.2f\n", *nparts, ComputeCut(&graph, part), lbvec[0], lbvec[1]);


  /*---------------------------------------------------------------------
   * Induce the decission tree
   *---------------------------------------------------------------------*/
  dtpart = idxmalloc(*nvtxs, "METIS_PartGraphForContact: dtpart");
  marker = idxsmalloc(*nvtxs, 0, "METIS_PartGraphForContact: marker");

  for (dim=0; dim<3; dim++) {
    xyzcand[dim] = (DKeyValueType *)gk_malloc(sizeof(DKeyValueType)*(*nvtxs), "METIS_PartGraphForContact: xyzcand[dim]");
    for (i=0; i<*nvtxs; i++) {
      xyzcand[dim][i].key = xyzcoords[3*i+dim];
      xyzcand[dim][i].val = i;
    }
    idkeysort(*nvtxs, xyzcand[dim]);
  }


  nnodes = nlnodes = nclean = naclean = ndirty = maxdepth = 0;
  InduceDecissionTree(*nvtxs, xyzcand, sflag, *nparts, part,
                      *nvtxs/(20*(*nparts)), *nvtxs/(20*(*nparts)*(*nparts)), 0.90,
                      &nnodes, &nlnodes, cinfo->dtree, leafpart, dtpart,
                      &nclean, &naclean, &ndirty, &maxdepth, marker);

  mprintf("NNodes: %5D, NLNodes: %5D, NClean: %5D, NAClean: %5D, NDirty: %5D, MaxDepth: %3D\n", nnodes, nlnodes, nclean, naclean, ndirty, maxdepth);


  /*---------------------------------------------------------------------
   * Create the tree-induced coarse graph and refine it
   *---------------------------------------------------------------------*/
  cgraph = CreatePartitionGraphForContact(*nvtxs, xadj, adjncy, mcvwgt, adjwgt, nlnodes, leafpart);

  for (i=0; i<*nvtxs; i++)
    part[leafpart[i]] = dtpart[i];

  ComputePartitionBalance(cgraph, *nparts, part, lbvec);
  mprintf("  %D-way Edge-Cut: %7D, Balance: %5.2f %5.2f\n", *nparts, ComputeCut(cgraph, part), lbvec[0], lbvec[1]);


  rwgtflag = 3;
  rnumflag = 0;
//.........这里部分代码省略.........

/*************************************************************************
* This function computes cuts and balance information
**************************************************************************/
void ComputePartitionInfo(GraphType *graph, int nparts, idxtype *where)
{
  int i, j, /*k,*/ nvtxs, ncon, mustfree=0;
  idxtype *xadj, *adjncy, *vwgt, *adjwgt, *kpwgts, *tmpptr;
  idxtype *padjncy, *padjwgt, *padjcut;

  nvtxs = graph->nvtxs;
  ncon = graph->ncon;
  xadj = graph->xadj;
  adjncy = graph->adjncy;
  vwgt = graph->vwgt;
  adjwgt = graph->adjwgt;

  if (vwgt == NULL) {
    vwgt = graph->vwgt = idxsmalloc(nvtxs, 1, "vwgt");
    mustfree = 1;
  }
  if (adjwgt == NULL) {
    adjwgt = graph->adjwgt = idxsmalloc(xadj[nvtxs], 1, "adjwgt");
    mustfree += 2;
  }

  printf("%d-way Cut: %5d, Vol: %5d, ", nparts, ComputeCut(graph, where), ComputeVolume(graph, where));

  /* Compute balance information */
  kpwgts = idxsmalloc(ncon*nparts, 0, "ComputePartitionInfo: kpwgts");

  for (i=0; i<nvtxs; i++) {
    for (j=0; j<ncon; j++) 
      kpwgts[where[i]*ncon+j] += vwgt[i*ncon+j];
  }

  if (ncon == 1) {
    printf("\tBalance: %5.3f out of %5.3f\n", 
            1.0*nparts*kpwgts[idxamax(nparts, kpwgts)]/(1.0*idxsum(nparts, kpwgts)),
            1.0*nparts*vwgt[idxamax(nvtxs, vwgt)]/(1.0*idxsum(nparts, kpwgts)));
  }
  else {
    printf("\tBalance:");
    for (j=0; j<ncon; j++) 
      printf(" (%5.3f out of %5.3f)", 
            1.0*nparts*kpwgts[ncon*idxamax_strd(nparts, kpwgts+j, ncon)+j]/(1.0*idxsum_strd(nparts, kpwgts+j, ncon)),
            1.0*nparts*vwgt[ncon*idxamax_strd(nvtxs, vwgt+j, ncon)+j]/(1.0*idxsum_strd(nparts, kpwgts+j, ncon)));
    printf("\n");
  }


  /* Compute p-adjncy information */
  padjncy = idxsmalloc(nparts*nparts, 0, "ComputePartitionInfo: padjncy");
  padjwgt = idxsmalloc(nparts*nparts, 0, "ComputePartitionInfo: padjwgt");
  padjcut = idxsmalloc(nparts*nparts, 0, "ComputePartitionInfo: padjwgt");

  idxset(nparts, 0, kpwgts);
  for (i=0; i<nvtxs; i++) {
    for (j=xadj[i]; j<xadj[i+1]; j++) {
      if (where[i] != where[adjncy[j]]) {
        padjncy[where[i]*nparts+where[adjncy[j]]] = 1;
        padjcut[where[i]*nparts+where[adjncy[j]]] += adjwgt[j];
        if (kpwgts[where[adjncy[j]]] == 0) {
          padjwgt[where[i]*nparts+where[adjncy[j]]]++;
          kpwgts[where[adjncy[j]]] = 1;
        }
      }
    }
    for (j=xadj[i]; j<xadj[i+1]; j++) 
      kpwgts[where[adjncy[j]]] = 0;
  }

  for (i=0; i<nparts; i++)
    kpwgts[i] = idxsum(nparts, padjncy+i*nparts);
  printf("Min/Max/Avg/Bal # of adjacent     subdomains: %5d %5d %5.2f %7.3f\n",
    kpwgts[idxamin(nparts, kpwgts)], kpwgts[idxamax(nparts, kpwgts)], 
    1.0*idxsum(nparts, kpwgts)/(1.0*nparts), 
    1.0*nparts*kpwgts[idxamax(nparts, kpwgts)]/(1.0*idxsum(nparts, kpwgts)));

  for (i=0; i<nparts; i++)
    kpwgts[i] = idxsum(nparts, padjcut+i*nparts);
  printf("Min/Max/Avg/Bal # of adjacent subdomain cuts: %5d %5d %5d %7.3f\n",
    kpwgts[idxamin(nparts, kpwgts)], kpwgts[idxamax(nparts, kpwgts)], idxsum(nparts, kpwgts)/nparts, 
    1.0*nparts*kpwgts[idxamax(nparts, kpwgts)]/(1.0*idxsum(nparts, kpwgts)));

  for (i=0; i<nparts; i++)
    kpwgts[i] = idxsum(nparts, padjwgt+i*nparts);
  printf("Min/Max/Avg/Bal/Frac # of interface    nodes: %5d %5d %5d %7.3f %7.3f\n",
    kpwgts[idxamin(nparts, kpwgts)], kpwgts[idxamax(nparts, kpwgts)], idxsum(nparts, kpwgts)/nparts, 
    1.0*nparts*kpwgts[idxamax(nparts, kpwgts)]/(1.0*idxsum(nparts, kpwgts)), 1.0*idxsum(nparts, kpwgts)/(1.0*nvtxs));

  tmpptr = graph->where;
  graph->where = where;
  for (i=0; i<nparts; i++)
    IsConnectedSubdomain(NULL, graph, i, 1);
  graph->where = tmpptr;

  if (mustfree == 1 || mustfree == 3) {
    free(vwgt);
    graph->vwgt = NULL;
  }
//.........这里部分代码省略.........

/*************************************************************************
* This function performs an edge-based FM refinement
**************************************************************************/
void MocGeneral2WayBalance(CtrlType *ctrl, GraphType *graph, float *tpwgts, float lbfactor)
{
  int i, ii, j, k, l, kwgt, nvtxs, ncon, nbnd, nswaps, from, to, pass, me, limit, tmp, cnum;
  idxtype *xadj, *adjncy, *adjwgt, *where, *id, *ed, *bndptr, *bndind;
  idxtype *moved, *swaps, *perm, *qnum;
  float *nvwgt, *npwgts, mindiff[MAXNCON], origbal, minbal, newbal;
  PQueueType parts[MAXNCON][2];
  int higain, oldgain, mincut, newcut, mincutorder;
  int qsizes[MAXNCON][2];

  nvtxs = graph->nvtxs;
  ncon = graph->ncon;
  xadj = graph->xadj;
  nvwgt = graph->nvwgt;
  adjncy = graph->adjncy;
  adjwgt = graph->adjwgt;
  where = graph->where;
  id = graph->id;
  ed = graph->ed;
  npwgts = graph->npwgts;
  bndptr = graph->bndptr;
  bndind = graph->bndind;

  moved = idxwspacemalloc(ctrl, nvtxs);
  swaps = idxwspacemalloc(ctrl, nvtxs);
  perm = idxwspacemalloc(ctrl, nvtxs);
  qnum = idxwspacemalloc(ctrl, nvtxs);

  limit = amin(amax(0.01*nvtxs, 15), 100);

  /* Initialize the queues */
  for (i=0; i<ncon; i++) {
    PQueueInit(ctrl, &parts[i][0], nvtxs, PLUS_GAINSPAN+1);
    PQueueInit(ctrl, &parts[i][1], nvtxs, PLUS_GAINSPAN+1);
    qsizes[i][0] = qsizes[i][1] = 0;
  }

  for (i=0; i<nvtxs; i++) {
    qnum[i] = samax(ncon, nvwgt+i*ncon);
    qsizes[qnum[i]][where[i]]++;
  }

/*
  printf("Weight Distribution:    \t");
  for (i=0; i<ncon; i++) 
    printf(" [%d %d]", qsizes[i][0], qsizes[i][1]); 
  printf("\n");
*/

  for (from=0; from<2; from++) {
    for (j=0; j<ncon; j++) {
      if (qsizes[j][from] == 0) {
        for (i=0; i<nvtxs; i++) {
          if (where[i] != from)
            continue;

          k = samax2(ncon, nvwgt+i*ncon);
          if (k == j && qsizes[qnum[i]][from] > qsizes[j][from] && nvwgt[i*ncon+qnum[i]] < 1.3*nvwgt[i*ncon+j]) {
            qsizes[qnum[i]][from]--;
            qsizes[j][from]++;
            qnum[i] = j;
          }
        }
      }
    }
  }

/*
  printf("Weight Distribution (after):\t ");
  for (i=0; i<ncon; i++) 
    printf(" [%d %d]", qsizes[i][0], qsizes[i][1]); 
  printf("\n");
*/



  for (i=0; i<ncon; i++) 
    mindiff[i] = fabs(tpwgts[0]-npwgts[i]);
  minbal = origbal = Compute2WayHLoadImbalance(ncon, npwgts, tpwgts);
  newcut = mincut = graph->mincut;
  mincutorder = -1;

  if (ctrl->dbglvl&DBG_REFINE) {
    printf("Parts: [");
    for (l=0; l<ncon; l++)
      printf("(%.3f, %.3f) ", npwgts[l], npwgts[ncon+l]);
    printf("] T[%.3f %.3f], Nv-Nb[%5d, %5d]. ICut: %6d, LB: %.3f [B]\n", tpwgts[0], tpwgts[1], graph->nvtxs, graph->nbnd, graph->mincut, origbal);
  }

  idxset(nvtxs, -1, moved);

  ASSERT(ComputeCut(graph, where) == graph->mincut);
  ASSERT(CheckBnd(graph));

  /* Insert all nodes in the priority queues */
  nbnd = graph->nbnd;
  RandomPermute(nvtxs, perm, 1);
//.........这里部分代码省略.........

/*************************************************************************
* This function performs an edge-based FM refinement
**************************************************************************/
void MocFM_2WayEdgeRefine(CtrlType *ctrl, GraphType *graph, float *tpwgts, int npasses)
{
  int i, ii, j, k, l, kwgt, nvtxs, ncon, nbnd, nswaps, from, to, pass, me, limit, tmp, cnum;
  idxtype *xadj, *adjncy, *adjwgt, *where, *id, *ed, *bndptr, *bndind;
  idxtype *moved, *swaps, *perm, *qnum;
  float *nvwgt, *npwgts, mindiff[MAXNCON], origbal, minbal, newbal;
  PQueueType parts[MAXNCON][2];
  int higain, oldgain, mincut, initcut, newcut, mincutorder;
  float rtpwgts[2];

  nvtxs = graph->nvtxs;
  ncon = graph->ncon;
  xadj = graph->xadj;
  nvwgt = graph->nvwgt;
  adjncy = graph->adjncy;
  adjwgt = graph->adjwgt;
  where = graph->where;
  id = graph->id;
  ed = graph->ed;
  npwgts = graph->npwgts;
  bndptr = graph->bndptr;
  bndind = graph->bndind;

  moved = idxwspacemalloc(ctrl, nvtxs);
  swaps = idxwspacemalloc(ctrl, nvtxs);
  perm = idxwspacemalloc(ctrl, nvtxs);
  qnum = idxwspacemalloc(ctrl, nvtxs);

  limit = amin(amax(0.01*nvtxs, 25), 150);

  /* Initialize the queues */
  for (i=0; i<ncon; i++) {
    PQueueInit(ctrl, &parts[i][0], nvtxs, PLUS_GAINSPAN+1);
    PQueueInit(ctrl, &parts[i][1], nvtxs, PLUS_GAINSPAN+1);
  }
  for (i=0; i<nvtxs; i++)
    qnum[i] = samax(ncon, nvwgt+i*ncon);

  origbal = Compute2WayHLoadImbalance(ncon, npwgts, tpwgts);

  rtpwgts[0] = origbal*tpwgts[0];
  rtpwgts[1] = origbal*tpwgts[1];

/*
  if (ctrl->dbglvl&DBG_REFINE) {
    printf("Parts: [");
    for (l=0; l<ncon; l++)
      printf("(%.3f, %.3f) ", npwgts[l], npwgts[ncon+l]);
    printf("] T[%.3f %.3f], Nv-Nb[%5d, %5d]. ICut: %6d, LB: %.3f\n", tpwgts[0], tpwgts[1], graph->nvtxs, graph->nbnd, graph->mincut, origbal);
  }
*/

  idxset(nvtxs, -1, moved);
  for (pass=0; pass<npasses; pass++) { /* Do a number of passes */
    for (i=0; i<ncon; i++) { 
      PQueueReset(&parts[i][0]);
      PQueueReset(&parts[i][1]);
    }

    mincutorder = -1;
    newcut = mincut = initcut = graph->mincut;
    for (i=0; i<ncon; i++)
      mindiff[i] = fabs(tpwgts[0]-npwgts[i]);
    minbal = Compute2WayHLoadImbalance(ncon, npwgts, tpwgts);

    ASSERT(ComputeCut(graph, where) == graph->mincut);
    ASSERT(CheckBnd(graph));

    /* Insert boundary nodes in the priority queues */
    nbnd = graph->nbnd;
    RandomPermute(nbnd, perm, 1);
    for (ii=0; ii<nbnd; ii++) {
      i = bndind[perm[ii]];
      ASSERT(ed[i] > 0 || id[i] == 0);
      ASSERT(bndptr[i] != -1);
      PQueueInsert(&parts[qnum[i]][where[i]], i, ed[i]-id[i]);
    }

    for (nswaps=0; nswaps<nvtxs; nswaps++) {
      SelectQueue(ncon, npwgts, rtpwgts, &from, &cnum, parts);
      to = (from+1)%2;

      if (from == -1 || (higain = PQueueGetMax(&parts[cnum][from])) == -1)
        break;
      ASSERT(bndptr[higain] != -1);

      saxpy(ncon, 1.0, nvwgt+higain*ncon, 1, npwgts+to*ncon, 1);
      saxpy(ncon, -1.0, nvwgt+higain*ncon, 1, npwgts+from*ncon, 1);

      newcut -= (ed[higain]-id[higain]);
      newbal = Compute2WayHLoadImbalance(ncon, npwgts, tpwgts);

      if ((newcut < mincut && newbal-origbal <= .00001) || 
          (newcut == mincut && (newbal < minbal || 
                                (newbal == minbal && BetterBalance(ncon, npwgts, tpwgts, mindiff))))) {
        mincut = newcut;
        minbal = newbal;
//.........这里部分代码省略.........

/*************************************************************************
* This function performs an edge-based FM refinement
**************************************************************************/
void FM_2WayEdgeRefine(CtrlType *ctrl, GraphType *graph, int *tpwgts, int npasses)
{
  int i, ii, j, k, kwgt, nvtxs, nbnd, nswaps, from, to, pass, me, limit, tmp;
  idxtype *xadj, *vwgt, *adjncy, *adjwgt, *where, *id, *ed, *bndptr, *bndind, *pwgts;
  idxtype *moved, *swaps, *perm;
  PQueueType parts[2];
  int higain, oldgain, mincut, mindiff, origdiff, initcut, newcut, mincutorder, avgvwgt;

  nvtxs = graph->nvtxs;
  xadj = graph->xadj;
  vwgt = graph->vwgt;
  adjncy = graph->adjncy;
  adjwgt = graph->adjwgt;
  where = graph->where;
  id = graph->id;
  ed = graph->ed;
  pwgts = graph->pwgts;
  bndptr = graph->bndptr;
  bndind = graph->bndind;

  moved = idxwspacemalloc(ctrl, nvtxs);
  swaps = idxwspacemalloc(ctrl, nvtxs);
  perm = idxwspacemalloc(ctrl, nvtxs);

  limit = (int) amin(amax(0.01*nvtxs, 15), 100);
  avgvwgt = amin((pwgts[0]+pwgts[1])/20, 2*(pwgts[0]+pwgts[1])/nvtxs);

  tmp = graph->adjwgtsum[idxamax(nvtxs, graph->adjwgtsum)];
  PQueueInit(ctrl, &parts[0], nvtxs, tmp);
  PQueueInit(ctrl, &parts[1], nvtxs, tmp);

  IFSET(ctrl->dbglvl, DBG_REFINE, 
     printf("Partitions: [%6d %6d] T[%6d %6d], Nv-Nb[%6d %6d]. ICut: %6d\n",
             pwgts[0], pwgts[1], tpwgts[0], tpwgts[1], graph->nvtxs, graph->nbnd, graph->mincut));

  origdiff = abs(tpwgts[0]-pwgts[0]);
  idxset(nvtxs, -1, moved);
  for (pass=0; pass<npasses; pass++) { /* Do a number of passes */
    PQueueReset(&parts[0]);
    PQueueReset(&parts[1]);

    mincutorder = -1;
    newcut = mincut = initcut = graph->mincut;
    mindiff = abs(tpwgts[0]-pwgts[0]);

    ASSERT(ComputeCut(graph, where) == graph->mincut);
    ASSERT(CheckBnd(graph));

    /* Insert boundary nodes in the priority queues */
    nbnd = graph->nbnd;
    RandomPermute(nbnd, perm, 1);
    for (ii=0; ii<nbnd; ii++) {
      i = perm[ii];
      ASSERT(ed[bndind[i]] > 0 || id[bndind[i]] == 0);
      ASSERT(bndptr[bndind[i]] != -1);
      PQueueInsert(&parts[where[bndind[i]]], bndind[i], ed[bndind[i]]-id[bndind[i]]);
    }

    for (nswaps=0; nswaps<nvtxs; nswaps++) {
      from = (tpwgts[0]-pwgts[0] < tpwgts[1]-pwgts[1] ? 0 : 1);
      to = (from+1)%2;

      if ((higain = PQueueGetMax(&parts[from])) == -1)
        break;
      ASSERT(bndptr[higain] != -1);

      newcut -= (ed[higain]-id[higain]);
      INC_DEC(pwgts[to], pwgts[from], vwgt[higain]);

      if ((newcut < mincut && abs(tpwgts[0]-pwgts[0]) <= origdiff+avgvwgt) || 
          (newcut == mincut && abs(tpwgts[0]-pwgts[0]) < mindiff)) {
        mincut = newcut;
        mindiff = abs(tpwgts[0]-pwgts[0]);
        mincutorder = nswaps;
      }
      else if (nswaps-mincutorder > limit) { /* We hit the limit, undo last move */
        newcut += (ed[higain]-id[higain]);
        INC_DEC(pwgts[from], pwgts[to], vwgt[higain]);
        break;
      }

      where[higain] = to;
      moved[higain] = nswaps;
      swaps[nswaps] = higain;

      IFSET(ctrl->dbglvl, DBG_MOVEINFO, 
        printf("Moved %6d from %d. [%3d %3d] %5d [%4d %4d]\n", higain, from, ed[higain]-id[higain], vwgt[higain], newcut, pwgts[0], pwgts[1]));

      /**************************************************************
      * Update the id[i]/ed[i] values of the affected nodes
      ***************************************************************/
      SWAP(id[higain], ed[higain], tmp);
      if (ed[higain] == 0 && xadj[higain] < xadj[higain+1]) 
        BNDDelete(nbnd, bndind,  bndptr, higain);

      for (j=xadj[higain]; j<xadj[higain+1]; j++) {
        k = adjncy[j];
//.........这里部分代码省略.........

/*************************************************************************
* This function balances two partitions by moving boundary nodes
* from the domain that is overweight to the one that is underweight.
**************************************************************************/
void Bnd2WayBalance(CtrlType *ctrl, GraphType *graph, int *tpwgts)
{
  int i, ii, j, k, kwgt, nvtxs, nbnd, nswaps, from, to, pass, me, tmp;
  idxtype *xadj, *vwgt, *adjncy, *adjwgt, *where, *id, *ed, *bndptr, *bndind, *pwgts;
  idxtype *moved, *perm;
  PQueueType parts;
  int higain, oldgain, mincut, mindiff;

  nvtxs = graph->nvtxs;
  xadj = graph->xadj;
  vwgt = graph->vwgt;
  adjncy = graph->adjncy;
  adjwgt = graph->adjwgt;
  where = graph->where;
  id = graph->id;
  ed = graph->ed;
  pwgts = graph->pwgts;
  bndptr = graph->bndptr;
  bndind = graph->bndind;

  moved = idxwspacemalloc(ctrl, nvtxs);
  perm = idxwspacemalloc(ctrl, nvtxs);

  /* Determine from which domain you will be moving data */
  mindiff = abs(tpwgts[0]-pwgts[0]);
  from = (pwgts[0] < tpwgts[0] ? 1 : 0);
  to = (from+1)%2;

  IFSET(ctrl->dbglvl, DBG_REFINE, 
     printf("Partitions: [%6d %6d] T[%6d %6d], Nv-Nb[%6d %6d]. ICut: %6d [B]\n",
             pwgts[0], pwgts[1], tpwgts[0], tpwgts[1], graph->nvtxs, graph->nbnd, graph->mincut));

  tmp = graph->adjwgtsum[idxamax(nvtxs, graph->adjwgtsum)];
  PQueueInit(ctrl, &parts, nvtxs, tmp);

  idxset(nvtxs, -1, moved);

  ASSERT(ComputeCut(graph, where) == graph->mincut);
  ASSERT(CheckBnd(graph));

  /* Insert the boundary nodes of the proper partition whose size is OK in the priority queue */
  nbnd = graph->nbnd;
  RandomPermute(nbnd, perm, 1);
  for (ii=0; ii<nbnd; ii++) {
    i = perm[ii];
    ASSERT(ed[bndind[i]] > 0 || id[bndind[i]] == 0);
    ASSERT(bndptr[bndind[i]] != -1);
    if (where[bndind[i]] == from && vwgt[bndind[i]] <= mindiff)
      PQueueInsert(&parts, bndind[i], ed[bndind[i]]-id[bndind[i]]);
  }

  mincut = graph->mincut;
  for (nswaps=0; nswaps<nvtxs; nswaps++) {
    if ((higain = PQueueGetMax(&parts)) == -1)
      break;
    ASSERT(bndptr[higain] != -1);

    if (pwgts[to]+vwgt[higain] > tpwgts[to])
      break;

    mincut -= (ed[higain]-id[higain]);
    INC_DEC(pwgts[to], pwgts[from], vwgt[higain]);

    where[higain] = to;
    moved[higain] = nswaps;

    IFSET(ctrl->dbglvl, DBG_MOVEINFO, 
      printf("Moved %6d from %d. [%3d %3d] %5d [%4d %4d]\n", higain, from, ed[higain]-id[higain], vwgt[higain], mincut, pwgts[0], pwgts[1]));

    /**************************************************************
    * Update the id[i]/ed[i] values of the affected nodes
    ***************************************************************/
    SWAP(id[higain], ed[higain], tmp);
    if (ed[higain] == 0 && xadj[higain] < xadj[higain+1]) 
      BNDDelete(nbnd, bndind,  bndptr, higain);

    for (j=xadj[higain]; j<xadj[higain+1]; j++) {
      k = adjncy[j];
      oldgain = ed[k]-id[k];

      kwgt = (to == where[k] ? adjwgt[j] : -adjwgt[j]);
      INC_DEC(id[k], ed[k], kwgt);

      /* Update its boundary information and queue position */
      if (bndptr[k] != -1) { /* If k was a boundary vertex */
        if (ed[k] == 0) { /* Not a boundary vertex any more */
          BNDDelete(nbnd, bndind, bndptr, k);
          if (moved[k] == -1 && where[k] == from && vwgt[k] <= mindiff)  /* Remove it if in the queues */
            PQueueDelete(&parts, k, oldgain);
        }
        else { /* If it has not been moved, update its position in the queue */
          if (moved[k] == -1 && where[k] == from && vwgt[k] <= mindiff)
            PQueueUpdate(&parts, k, oldgain, ed[k]-id[k]);
        }
      }
      else {
//.........这里部分代码省略.........

                        }
                    }
                }
            }
            graph->minvol -= xgain;
            graph->mincut -= -myrinfo->nid;
        } else {
            graph->minvol -= (xgain + mynbrs[k].gv);
            graph->mincut -= mynbrs[k].ned - myrinfo->nid;
        }


        /* Update where and pwgts */
        where[i] = to;
        iaxpy(graph->ncon, 1, graph->vwgt + i * graph->ncon, 1, graph->pwgts + to * graph->ncon, 1);
        iaxpy(graph->ncon, -1, graph->vwgt + i * graph->ncon, 1, graph->pwgts + from * graph->ncon, 1);

        /* Update the id/ed/gains/bnd of potentially affected nodes */
        KWayVolUpdate(ctrl, graph, i, from, to, NULL, NULL, NULL, NULL,
                      NULL, BNDTYPE_REFINE, vmarker, pmarker, modind);

        /*CheckKWayVolPartitionParams(ctrl, graph);*/
    }

    ASSERT(ComputeCut(graph, where) == graph->mincut);
    ASSERTP(ComputeVolume(graph, where) == graph->minvol,
            ("%"PRIDX" %"PRIDX"\n", ComputeVolume(graph, where), graph->minvol));

}

/*************************************************************************
* This function computes cuts and balance information
**************************************************************************/
void ComputePartitionInfoBipartite(graph_t *graph, idx_t nparts, idx_t *where)
{
  idx_t i, j, k, nvtxs, ncon, mustfree=0;
  idx_t *xadj, *adjncy, *vwgt, *vsize, *adjwgt, *kpwgts, *tmpptr;
  idx_t *padjncy, *padjwgt, *padjcut;

  nvtxs = graph->nvtxs;
  ncon = graph->ncon;
  xadj = graph->xadj;
  adjncy = graph->adjncy;
  vwgt = graph->vwgt;
  vsize = graph->vsize;
  adjwgt = graph->adjwgt;

  if (vwgt == NULL) {
    vwgt = graph->vwgt = ismalloc(nvtxs, 1, "vwgt");
    mustfree = 1;
  }
  if (adjwgt == NULL) {
    adjwgt = graph->adjwgt = ismalloc(xadj[nvtxs], 1, "adjwgt");
    mustfree += 2;
  }

  printf("%"PRIDX"-way Cut: %5"PRIDX", Vol: %5"PRIDX", ", nparts, ComputeCut(graph, where), ComputeVolume(graph, where));

  /* Compute balance information */
  kpwgts = ismalloc(ncon*nparts, 0, "ComputePartitionInfo: kpwgts");

  for (i=0; i<nvtxs; i++) {
    for (j=0; j<ncon; j++) 
      kpwgts[where[i]*ncon+j] += vwgt[i*ncon+j];
  }

  if (ncon == 1) {
    printf("\tBalance: %5.3"PRREAL" out of %5.3"PRREAL"\n", 
            1.0*nparts*kpwgts[iargmax(nparts, kpwgts)]/(1.0*isum(nparts, kpwgts, 1)),
            1.0*nparts*vwgt[iargmax(nvtxs, vwgt)]/(1.0*isum(nparts, kpwgts, 1)));
  }
  else {
    printf("\tBalance:");
    for (j=0; j<ncon; j++) 
      printf(" (%5.3"PRREAL" out of %5.3"PRREAL")", 
            1.0*nparts*kpwgts[ncon*iargmax_strd(nparts, kpwgts+j, ncon)+j]/(1.0*isum(nparts, kpwgts+j, ncon)),
            1.0*nparts*vwgt[ncon*iargmax_strd(nvtxs, vwgt+j, ncon)+j]/(1.0*isum(nparts, kpwgts+j, ncon)));
    printf("\n");
  }


  /* Compute p-adjncy information */
  padjncy = ismalloc(nparts*nparts, 0, "ComputePartitionInfo: padjncy");
  padjwgt = ismalloc(nparts*nparts, 0, "ComputePartitionInfo: padjwgt");
  padjcut = ismalloc(nparts*nparts, 0, "ComputePartitionInfo: padjwgt");

  iset(nparts, 0, kpwgts);
  for (i=0; i<nvtxs; i++) {
    for (j=xadj[i]; j<xadj[i+1]; j++) {
      if (where[i] != where[adjncy[j]]) {
        padjncy[where[i]*nparts+where[adjncy[j]]] = 1;
        padjcut[where[i]*nparts+where[adjncy[j]]] += adjwgt[j];
        if (kpwgts[where[adjncy[j]]] == 0) {
          padjwgt[where[i]*nparts+where[adjncy[j]]] += vsize[i];
          kpwgts[where[adjncy[j]]] = 1;
        }
      }
    }
    for (j=xadj[i]; j<xadj[i+1]; j++) 
      kpwgts[where[adjncy[j]]] = 0;
  }

  for (i=0; i<nparts; i++)
    kpwgts[i] = isum(nparts, padjncy+i*nparts, 1);
  printf("Min/Max/Avg/Bal # of adjacent     subdomains: %5"PRIDX" %5"PRIDX" %5"PRIDX" %7.3"PRREAL"\n",
    kpwgts[iargmin(nparts, kpwgts)], kpwgts[iargmax(nparts, kpwgts)], isum(nparts, kpwgts, 1)/nparts, 
    1.0*nparts*kpwgts[iargmax(nparts, kpwgts)]/(1.0*isum(nparts, kpwgts, 1)));

  for (i=0; i<nparts; i++)
    kpwgts[i] = isum(nparts, padjcut+i*nparts, 1);
  printf("Min/Max/Avg/Bal # of adjacent subdomain cuts: %5"PRIDX" %5"PRIDX" %5"PRIDX" %7.3"PRREAL"\n",
    kpwgts[iargmin(nparts, kpwgts)], kpwgts[iargmax(nparts, kpwgts)], isum(nparts, kpwgts, 1)/nparts, 
    1.0*nparts*kpwgts[iargmax(nparts, kpwgts)]/(1.0*isum(nparts, kpwgts, 1)));

  for (i=0; i<nparts; i++)
    kpwgts[i] = isum(nparts, padjwgt+i*nparts, 1);
  printf("Min/Max/Avg/Bal/Frac # of interface    nodes: %5"PRIDX" %5"PRIDX" %5"PRIDX" %7.3"PRREAL" %7.3"PRREAL"\n",
    kpwgts[iargmin(nparts, kpwgts)], kpwgts[iargmax(nparts, kpwgts)], isum(nparts, kpwgts, 1)/nparts, 
    1.0*nparts*kpwgts[iargmax(nparts, kpwgts)]/(1.0*isum(nparts, kpwgts, 1)), 1.0*isum(nparts, kpwgts, 1)/(1.0*nvtxs));


  if (mustfree == 1 || mustfree == 3) {
    gk_free((void **)&vwgt, LTERM);
    graph->vwgt = NULL;
  }
  if (mustfree == 2 || mustfree == 3) {
    gk_free((void **)&adjwgt, LTERM);
    graph->adjwgt = NULL;
  }

//.........这里部分代码省略.........

/*************************************************************************
* This function performs k-way refinement
**************************************************************************/
void Random_KWayEdgeRefineMConn(CtrlType *ctrl, GraphType *graph, int nparts, float *tpwgts, float ubfactor, int npasses, int ffactor)
{
  int i, ii, iii, j, jj, k, l, pass, nvtxs, nmoves, nbnd, tvwgt, myndegrees; 
  int from, me, to, oldcut, vwgt, gain;
  int maxndoms, nadd;
  idxtype *xadj, *adjncy, *adjwgt;
  idxtype *where, *pwgts, *perm, *bndptr, *bndind, *minwgt, *maxwgt, *itpwgts;
  idxtype *phtable, *pmat, *pmatptr, *ndoms;
  EDegreeType *myedegrees;
  RInfoType *myrinfo;

  nvtxs = graph->nvtxs;
  xadj = graph->xadj;
  adjncy = graph->adjncy;
  adjwgt = graph->adjwgt;

  bndptr = graph->bndptr;
  bndind = graph->bndind;

  where = graph->where;
  pwgts = graph->pwgts;

  pmat = ctrl->wspace.pmat;
  phtable = idxwspacemalloc(ctrl, nparts);
  ndoms = idxwspacemalloc(ctrl, nparts);

  ComputeSubDomainGraph(graph, nparts, pmat, ndoms);

  /* Setup the weight intervals of the various subdomains */
  minwgt =  idxwspacemalloc(ctrl, nparts);
  maxwgt = idxwspacemalloc(ctrl, nparts);
  itpwgts = idxwspacemalloc(ctrl, nparts);
  tvwgt = idxsum(nparts, pwgts);
  ASSERT(tvwgt == idxsum(nvtxs, graph->vwgt));

  for (i=0; i<nparts; i++) {
    itpwgts[i] = tpwgts[i]*tvwgt;
    maxwgt[i] = tpwgts[i]*tvwgt*ubfactor;
    minwgt[i] = tpwgts[i]*tvwgt*(1.0/ubfactor);
  }

  perm = idxwspacemalloc(ctrl, nvtxs);

  IFSET(ctrl->dbglvl, DBG_REFINE,
     printf("Partitions: [%6d %6d]-[%6d %6d], Balance: %5.3f, Nv-Nb[%6d %6d]. Cut: %6d\n",
             pwgts[idxamin(nparts, pwgts)], pwgts[idxamax(nparts, pwgts)], minwgt[0], maxwgt[0], 
             1.0*nparts*pwgts[idxamax(nparts, pwgts)]/tvwgt, graph->nvtxs, graph->nbnd,
             graph->mincut));

  for (pass=0; pass<npasses; pass++) {
    ASSERT(ComputeCut(graph, where) == graph->mincut);

    maxndoms = ndoms[idxamax(nparts, ndoms)];

    oldcut = graph->mincut;
    nbnd = graph->nbnd;

    RandomPermute(nbnd, perm, 1);
    for (nmoves=iii=0; iii<graph->nbnd; iii++) {
      ii = perm[iii];
      if (ii >= nbnd)
        continue;
      i = bndind[ii];

      myrinfo = graph->rinfo+i;

      if (myrinfo->ed >= myrinfo->id) { /* Total ED is too high */
        from = where[i];
        vwgt = graph->vwgt[i];

        if (myrinfo->id > 0 && pwgts[from]-vwgt < minwgt[from]) 
          continue;   /* This cannot be moved! */

        myedegrees = myrinfo->edegrees;
        myndegrees = myrinfo->ndegrees;

        /* Determine the valid domains */
        for (j=0; j<myndegrees; j++) {
          to = myedegrees[j].pid;
          phtable[to] = 1;
          pmatptr = pmat + to*nparts;
          for (nadd=0, k=0; k<myndegrees; k++) {
            if (k == j)
              continue;

            l = myedegrees[k].pid;
            if (pmatptr[l] == 0) {
              if (ndoms[l] > maxndoms-1) {
                phtable[to] = 0;
                nadd = maxndoms;
                break;
              }
              nadd++;
            }
          }
          if (ndoms[to]+nadd > maxndoms)
            phtable[to] = 0;
//.........这里部分代码省略.........

/*************************************************************************
* This function performs k-way refinement
**************************************************************************/
void Greedy_KWayEdgeRefine(CtrlType *ctrl, GraphType *graph, int nparts, float *tpwgts, float ubfactor, int npasses)
{
    int i, ii, iii, j, jj, k, l, pass, nvtxs, nbnd, tvwgt, myndegrees, oldgain, gain;
    int from, me, to, oldcut, vwgt;
    idxtype *xadj, *adjncy, *adjwgt;
    idxtype *where, *pwgts, *perm, *bndptr, *bndind, *minwgt, *maxwgt, *moved, *itpwgts;
    EDegreeType *myedegrees;
    RInfoType *myrinfo;
    PQueueType queue;

    nvtxs = graph->nvtxs;
    xadj = graph->xadj;
    adjncy = graph->adjncy;
    adjwgt = graph->adjwgt;

    bndind = graph->bndind;
    bndptr = graph->bndptr;

    where = graph->where;
    pwgts = graph->pwgts;

    /* Setup the weight intervals of the various subdomains */
    minwgt =  idxwspacemalloc(ctrl, nparts);
    maxwgt = idxwspacemalloc(ctrl, nparts);
    itpwgts = idxwspacemalloc(ctrl, nparts);
    tvwgt = idxsum(nparts, pwgts);
    ASSERT(tvwgt == idxsum(nvtxs, graph->vwgt));

    for (i=0; i<nparts; i++) {
        itpwgts[i] = tpwgts[i]*tvwgt;
        maxwgt[i] = tpwgts[i]*tvwgt*ubfactor;
        minwgt[i] = tpwgts[i]*tvwgt*(1.0/ubfactor);
    }

    perm = idxwspacemalloc(ctrl, nvtxs);
    moved = idxwspacemalloc(ctrl, nvtxs);

    PQueueInit(ctrl, &queue, nvtxs, graph->adjwgtsum[idxamax(nvtxs, graph->adjwgtsum)]);

    IFSET(ctrl->dbglvl, DBG_REFINE,
          printf("Partitions: [%6d %6d]-[%6d %6d], Balance: %5.3f, Nv-Nb[%6d %6d]. Cut: %6d\n",
                 pwgts[idxamin(nparts, pwgts)], pwgts[idxamax(nparts, pwgts)], minwgt[0], maxwgt[0],
                 1.0*nparts*pwgts[idxamax(nparts, pwgts)]/tvwgt, graph->nvtxs, graph->nbnd,
                 graph->mincut));

    for (pass=0; pass<npasses; pass++) {
        ASSERT(ComputeCut(graph, where) == graph->mincut);

        PQueueReset(&queue);
        idxset(nvtxs, -1, moved);

        oldcut = graph->mincut;
        nbnd = graph->nbnd;

        RandomPermute(nbnd, perm, 1);
        for (ii=0; ii<nbnd; ii++) {
            i = bndind[perm[ii]];
            PQueueInsert(&queue, i, graph->rinfo[i].ed - graph->rinfo[i].id);
            moved[i] = 2;
        }

        for (iii=0;; iii++) {
            if ((i = PQueueGetMax(&queue)) == -1)
                break;
            moved[i] = 1;

            myrinfo = graph->rinfo+i;
            from = where[i];
            vwgt = graph->vwgt[i];

            if (pwgts[from]-vwgt < minwgt[from])
                continue;   /* This cannot be moved! */

            myedegrees = myrinfo->edegrees;
            myndegrees = myrinfo->ndegrees;

            j = myrinfo->id;
            for (k=0; k<myndegrees; k++) {
                to = myedegrees[k].pid;
                gain = myedegrees[k].ed-j; /* j = myrinfo->id. Allow good nodes to move */
                if (pwgts[to]+vwgt <= maxwgt[to]+gain && gain >= 0)
                    break;
            }
            if (k == myndegrees)
                continue;  /* break out if you did not find a candidate */

            for (j=k+1; j<myndegrees; j++) {
                to = myedegrees[j].pid;
                if ((myedegrees[j].ed > myedegrees[k].ed && pwgts[to]+vwgt <= maxwgt[to]) ||
                        (myedegrees[j].ed == myedegrees[k].ed &&
                         itpwgts[myedegrees[k].pid]*pwgts[to] < itpwgts[to]*pwgts[myedegrees[k].pid]))
                    k = j;
            }

            to = myedegrees[k].pid;

            j = 0;
//.........这里部分代码省略.........