src/libs/ck-libs/metis/libmetis/parmetis.c

   1 /*
   2  * Copyright 1997, Regents of the University of Minnesota
   3  *
   4  * parmetis.c
   5  *
   6  * This file contains top level routines that are used by ParMETIS
   7  *
   8  * Started 10/14/97
   9  * George
  10  *
  11  * $Id: parmetis.c 10481 2011-07-05 18:01:23Z karypis $
  12  *
  13  */
  14
  15 #include "metislib.h"
  16
  17
  18 /*************************************************************************/
  19 /*! This function is the entry point for the node ND code for ParMETIS.
  20     The difference between this routine and the standard METIS_NodeND are
  21     the following
  22
  23     - It performs at least log2(npes) levels of nested dissection.
  24     - It stores the size of the log2(npes) top-level separators in the
  25       sizes array.
  26 */
  27 /*************************************************************************/
  28 int METIS_NodeNDP(idx_t nvtxs, idx_t *xadj, idx_t *adjncy, idx_t *vwgt,
  29            idx_t npes, idx_t *options, idx_t *perm, idx_t *iperm, idx_t *sizes)
  30 {
  31   idx_t i, ii, j, l, nnvtxs=0;
  32   graph_t *graph;
  33   ctrl_t *ctrl;
  34   idx_t *cptr, *cind;
  35
  36   ctrl = SetupCtrl(METIS_OP_OMETIS, options, 1, 3, NULL, NULL);
  37   if (!ctrl) return METIS_ERROR_INPUT;
  38
  39   IFSET(ctrl->dbglvl, METIS_DBG_TIME, InitTimers(ctrl));
  40   IFSET(ctrl->dbglvl, METIS_DBG_TIME, gk_startcputimer(ctrl->TotalTmr));
  41
  42   /* compress the graph; not that compression only happens if not prunning
  43      has taken place. */
  44   if (ctrl->compress) {
  45     cptr = imalloc(nvtxs+1, "OMETIS: cptr");
  46     cind = imalloc(nvtxs, "OMETIS: cind");
  47
  48     graph = CompressGraph(ctrl, nvtxs, xadj, adjncy, vwgt, cptr, cind);
  49     if (graph == NULL) {
  50       /* if there was no compression, cleanup the compress flag */
  51       gk_free((void **)&cptr, &cind, LTERM);
  52       ctrl->compress = 0;
  53     }
  54     else {
  55       nnvtxs = graph->nvtxs;
  56     }
  57   }
  58
  59   /* if no compression, setup the graph in the normal way. */
  60   if (ctrl->compress == 0)
  61     graph = SetupGraph(ctrl, nvtxs, 1, xadj, adjncy, vwgt, NULL, NULL);
  62
  63
  64   /* allocate workspace memory */
  65   AllocateWorkSpace(ctrl, graph);
  66
  67
  68   /* do the nested dissection ordering  */
  69   iset(2*npes-1, 0, sizes);
  70   MlevelNestedDissectionP(ctrl, graph, iperm, graph->nvtxs, npes, 0, sizes);
  71
  72
  73   /* Uncompress the ordering */
  74   if (ctrl->compress) {
  75     /* construct perm from iperm */
  76     for (i=0; i<nnvtxs; i++)
  77       perm[iperm[i]] = i;
  78     for (l=ii=0; ii<nnvtxs; ii++) {
  79       i = perm[ii];
  80       for (j=cptr[i]; j<cptr[i+1]; j++)
  81         iperm[cind[j]] = l++;
  82     }
  83
  84     gk_free((void **)&cptr, &cind, LTERM);
  85   }
  86
  87
  88   for (i=0; i<nvtxs; i++)
  89     perm[iperm[i]] = i;
  90
  91   IFSET(ctrl->dbglvl, METIS_DBG_TIME, gk_stopcputimer(ctrl->TotalTmr));
  92   IFSET(ctrl->dbglvl, METIS_DBG_TIME, PrintTimers(ctrl));
  93
  94   /* clean up */
  95   FreeCtrl(&ctrl);
  96
  97   return METIS_OK;
  98 }
  99
 100
 101 /*************************************************************************/
 102 /*! This function is similar to MlevelNestedDissection with the difference
 103     that it also records separator sizes for the top log2(npes) levels */
 104 /**************************************************************************/
 105 void MlevelNestedDissectionP(ctrl_t *ctrl, graph_t *graph, idx_t *order,
 106          idx_t lastvtx, idx_t npes, idx_t cpos, idx_t *sizes)
 107 {
 108   idx_t i, j, nvtxs, nbnd;
 109   idx_t *label, *bndind;
 110   graph_t *lgraph, *rgraph;
 111
 112   nvtxs = graph->nvtxs;
 113
 114   if (nvtxs == 0) {
 115     FreeGraph(&graph);
 116     return;
 117   }
 118
 119   MlevelNodeBisectionMultiple(ctrl, graph);
 120
 121   IFSET(ctrl->dbglvl, METIS_DBG_SEPINFO,
 122       printf("Nvtxs: %6"PRIDX", [%6"PRIDX" %6"PRIDX" %6"PRIDX"]\n",
 123         graph->nvtxs, graph->pwgts[0], graph->pwgts[1], graph->pwgts[2]));
 124
 125   if (cpos < npes-1) {
 126     sizes[2*npes-2-cpos]       = graph->pwgts[2];
 127     sizes[2*npes-2-(2*cpos+1)] = graph->pwgts[1];
 128     sizes[2*npes-2-(2*cpos+2)] = graph->pwgts[0];
 129   }
 130
 131   /* Order the nodes in the separator */
 132   nbnd   = graph->nbnd;
 133   bndind = graph->bndind;
 134   label  = graph->label;
 135   for (i=0; i<nbnd; i++)
 136     order[label[bndind[i]]] = --lastvtx;
 137
 138   SplitGraphOrder(ctrl, graph, &lgraph, &rgraph);
 139
 140   /* Free the memory of the top level graph */
 141   FreeGraph(&graph);
 142
 143   if ((lgraph->nvtxs > MMDSWITCH || 2*cpos+2 < npes-1) && lgraph->nedges > 0)
 144     MlevelNestedDissectionP(ctrl, lgraph, order, lastvtx-rgraph->nvtxs, npes, 2*cpos+2, sizes);
 145   else {
 146     MMDOrder(ctrl, lgraph, order, lastvtx-rgraph->nvtxs);
 147     FreeGraph(&lgraph);
 148   }
 149   if ((rgraph->nvtxs > MMDSWITCH || 2*cpos+1 < npes-1) && rgraph->nedges > 0)
 150     MlevelNestedDissectionP(ctrl, rgraph, order, lastvtx, npes, 2*cpos+1, sizes);
 151   else {
 152     MMDOrder(ctrl, rgraph, order, lastvtx);
 153     FreeGraph(&rgraph);
 154   }
 155 }
 156
 157
 158 /*************************************************************************/
 159 /*! This function bisects a graph by computing a vertex separator */
 160 /**************************************************************************/
 161 int METIS_ComputeVertexSeparator(idx_t *nvtxs, idx_t *xadj, idx_t *adjncy,
 162            idx_t *vwgt, idx_t *options, idx_t *r_sepsize, idx_t *part)
 163 {
 164   idx_t i, j;
 165   graph_t *graph;
 166   ctrl_t *ctrl;
 167
 168   if ((ctrl = SetupCtrl(METIS_OP_OMETIS, options, 1, 3, NULL, NULL)) == NULL)
 169     return METIS_ERROR_INPUT;
 170
 171   InitRandom(ctrl->seed);
 172
 173   graph = SetupGraph(ctrl, *nvtxs, 1, xadj, adjncy, vwgt, NULL, NULL);
 174
 175   AllocateWorkSpace(ctrl, graph);
 176
 177   /*============================================================
 178    * Perform the bisection
 179    *============================================================*/
 180   ctrl->CoarsenTo = 100;
 181
 182   MlevelNodeBisectionMultiple(ctrl, graph);
 183
 184   *r_sepsize = graph->pwgts[2];
 185   icopy(*nvtxs, graph->where, part);
 186
 187   FreeGraph(&graph);
 188
 189   FreeCtrl(&ctrl);
 190
 191   return METIS_OK;
 192 }
 193
 194
 195 /*************************************************************************/
 196 /*! This function is the entry point of a node-based separator refinement
 197     of the nodes with an hmarker[] of 0. */
 198 /*************************************************************************/
 199 int METIS_NodeRefine(idx_t nvtxs, idx_t *xadj, idx_t *vwgt, idx_t *adjncy,
 200            idx_t *where, idx_t *hmarker, real_t ubfactor)
 201 {
 202   graph_t *graph;
 203   ctrl_t *ctrl;
 204
 205   /* set up the run time parameters */
 206   ctrl = SetupCtrl(METIS_OP_OMETIS, NULL, 1, 3, NULL, NULL);
 207   if (!ctrl) return METIS_ERROR_INPUT;
 208
 209   /* set up the graph */
 210   graph = SetupGraph(ctrl, nvtxs, 1, xadj, adjncy, vwgt, NULL, NULL);
 211
 212   /* allocate workspace memory */
 213   AllocateWorkSpace(ctrl, graph);
 214
 215   /* set up the memory and the input partition */
 216   Allocate2WayNodePartitionMemory(ctrl, graph);
 217   icopy(nvtxs, where, graph->where);
 218
 219   Compute2WayNodePartitionParams(ctrl, graph);
 220
 221   FM_2WayNodeRefine1SidedP(ctrl, graph, hmarker, ubfactor, 10);
 222   /* FM_2WayNodeRefine2SidedP(ctrl, graph, hmarker, ubfactor, 10); */
 223
 224   icopy(nvtxs, graph->where, where);
 225
 226   FreeGraph(&graph);
 227   FreeCtrl(&ctrl);
 228
 229   return METIS_OK;
 230 }
 231
 232
 233 /*************************************************************************/
 234 /*! This function performs a node-based 1-sided FM refinement that moves
 235     only nodes whose hmarker[] == -1. It is used by Parmetis. */
 236 /*************************************************************************/
 237 void FM_2WayNodeRefine1SidedP(ctrl_t *ctrl, graph_t *graph,
 238           idx_t *hmarker, real_t ubfactor, idx_t npasses)
 239 {
 240   idx_t i, ii, j, k, jj, kk, nvtxs, nbnd, nswaps, nmind, nbad, qsize;
 241   idx_t *xadj, *vwgt, *adjncy, *where, *pwgts, *edegrees, *bndind, *bndptr;
 242   idx_t *mptr, *mind, *swaps, *inqueue;
 243   rpq_t *queue;
 244   nrinfo_t *rinfo;
 245   idx_t higain, oldgain, mincut, initcut, mincutorder;
 246   idx_t pass, from, to, limit;
 247   idx_t badmaxpwgt, mindiff, newdiff;
 248
 249   WCOREPUSH;
 250
 251   ASSERT(graph->mincut == graph->pwgts[2]);
 252
 253   nvtxs  = graph->nvtxs;
 254   xadj   = graph->xadj;
 255   adjncy = graph->adjncy;
 256   vwgt   = graph->vwgt;
 257
 258   bndind = graph->bndind;
 259   bndptr = graph->bndptr;
 260   where  = graph->where;
 261   pwgts  = graph->pwgts;
 262   rinfo  = graph->nrinfo;
 263
 264   queue = rpqCreate(nvtxs);
 265
 266   inqueue = iset(nvtxs, -1, iwspacemalloc(ctrl, nvtxs));
 267   swaps   = iwspacemalloc(ctrl, nvtxs);
 268   mptr    = iwspacemalloc(ctrl, nvtxs+1);
 269   mind    = iwspacemalloc(ctrl, 2*nvtxs);
 270
 271   badmaxpwgt = (idx_t)(ubfactor*gk_max(pwgts[0], pwgts[1]));
 272
 273   IFSET(ctrl->dbglvl, METIS_DBG_REFINE,
 274     printf("Partitions-N1: [%6"PRIDX" %6"PRIDX"] Nv-Nb[%6"PRIDX" %6"PRIDX"] "
 275            "MaxPwgt[%6"PRIDX"]. ISep: %6"PRIDX"\n",
 276            pwgts[0], pwgts[1], graph->nvtxs, graph->nbnd, badmaxpwgt,
 277            graph->mincut));
 278
 279   to = (pwgts[0] < pwgts[1] ? 1 : 0);
 280   for (pass=0; pass<npasses; pass++) {
 281     from = to;
 282     to   = (from+1)%2;
 283
 284     rpqReset(queue);
 285
 286     mincutorder = -1;
 287     initcut = mincut = graph->mincut;
 288     nbnd = graph->nbnd;
 289
 290     /* use the swaps array in place of the traditional perm array to save memory */
 291     irandArrayPermute(nbnd, swaps, nbnd, 1);
 292     for (ii=0; ii<nbnd; ii++) {
 293       i = bndind[swaps[ii]];
 294       ASSERT(where[i] == 2);
 295       if (hmarker[i] == -1 || hmarker[i] == to) {
 296         rpqInsert(queue, i, vwgt[i]-rinfo[i].edegrees[from]);
 297         inqueue[i] = pass;
 298       }
 299     }
 300     qsize = rpqLength(queue);
 301
 302     ASSERT(CheckNodeBnd(graph, nbnd));
 303     ASSERT(CheckNodePartitionParams(graph));
 304
 305     limit = nbnd;
 306
 307     /******************************************************
 308     * Get into the FM loop
 309     *******************************************************/
 310     mptr[0] = nmind = nbad = 0;
 311     mindiff = abs(pwgts[0]-pwgts[1]);
 312     for (nswaps=0; nswaps<nvtxs; nswaps++) {
 313       if ((higain = rpqGetTop(queue)) == -1)
 314         break;
 315
 316       ASSERT(bndptr[higain] != -1);
 317
 318       /* The following check is to ensure we break out if there is a posibility
 319          of over-running the mind array.  */
 320       if (nmind + xadj[higain+1]-xadj[higain] >= 2*nvtxs-1)
 321         break;
 322
 323       inqueue[higain] = -1;
 324
 325       if (pwgts[to]+vwgt[higain] > badmaxpwgt) { /* Skip this vertex */
 326         if (nbad++ > limit)
 327           break;
 328         else {
 329           nswaps--;
 330           continue;
 331         }
 332       }
 333
 334       pwgts[2] -= (vwgt[higain]-rinfo[higain].edegrees[from]);
 335
 336       newdiff = abs(pwgts[to]+vwgt[higain] - (pwgts[from]-rinfo[higain].edegrees[from]));
 337       if (pwgts[2] < mincut || (pwgts[2] == mincut && newdiff < mindiff)) {
 338         mincut      = pwgts[2];
 339         mincutorder = nswaps;
 340         mindiff     = newdiff;
 341         nbad        = 0;
 342       }
 343       else {
 344         if (nbad++ > limit) {
 345           pwgts[2] += (vwgt[higain]-rinfo[higain].edegrees[from]);
 346           break; /* No further improvement, break out */
 347         }
 348       }
 349
 350       BNDDelete(nbnd, bndind, bndptr, higain);
 351       pwgts[to] += vwgt[higain];
 352       where[higain] = to;
 353       swaps[nswaps] = higain;
 354
 355
 356       /**********************************************************
 357       * Update the degrees of the affected nodes
 358       ***********************************************************/
 359       for (j=xadj[higain]; j<xadj[higain+1]; j++) {
 360         k = adjncy[j];
 361         if (where[k] == 2) { /* For the in-separator vertices modify their edegree[to] */
 362           rinfo[k].edegrees[to] += vwgt[higain];
 363         }
 364         else if (where[k] == from) { /* This vertex is pulled into the separator */
 365           ASSERTP(bndptr[k] == -1, ("%"PRIDX" %"PRIDX" %"PRIDX"\n", k, bndptr[k], where[k]));
 366           BNDInsert(nbnd, bndind, bndptr, k);
 367
 368           mind[nmind++] = k;  /* Keep track for rollback */
 369           where[k]      = 2;
 370           pwgts[from]  -= vwgt[k];
 371
 372           edegrees = rinfo[k].edegrees;
 373           edegrees[0] = edegrees[1] = 0;
 374           for (jj=xadj[k]; jj<xadj[k+1]; jj++) {
 375             kk = adjncy[jj];
 376             if (where[kk] != 2)
 377               edegrees[where[kk]] += vwgt[kk];
 378             else {
 379               oldgain = vwgt[kk]-rinfo[kk].edegrees[from];
 380               rinfo[kk].edegrees[from] -= vwgt[k];
 381
 382               /* Update the gain of this node if it was not skipped */
 383               if (inqueue[kk] == pass)
 384                 rpqUpdate(queue, kk, oldgain+vwgt[k]);
 385             }
 386           }
 387
 388           /* Insert the new vertex into the priority queue. Safe due to one-sided moves */
 389           if (hmarker[k] == -1 || hmarker[k] == to) {
 390             rpqInsert(queue, k, vwgt[k]-edegrees[from]);
 391             inqueue[k] = pass;
 392           }
 393         }
 394       }
 395       mptr[nswaps+1] = nmind;
 396
 397
 398       IFSET(ctrl->dbglvl, METIS_DBG_MOVEINFO,
 399             printf("Moved %6"PRIDX" to %3"PRIDX", Gain: %5"PRIDX" [%5"PRIDX"] \t[%5"PRIDX" %5"PRIDX" %5"PRIDX"] [%3"PRIDX" %2"PRIDX"]\n",
 400                    higain, to, (vwgt[higain]-rinfo[higain].edegrees[from]),
 401                    vwgt[higain], pwgts[0], pwgts[1], pwgts[2], nswaps, limit));
 402
 403     }
 404
 405
 406     /****************************************************************
 407     * Roll back computation
 408     *****************************************************************/
 409     for (nswaps--; nswaps>mincutorder; nswaps--) {
 410       higain = swaps[nswaps];
 411
 412       ASSERT(CheckNodePartitionParams(graph));
 413       ASSERT(where[higain] == to);
 414
 415       INC_DEC(pwgts[2], pwgts[to], vwgt[higain]);
 416       where[higain] = 2;
 417       BNDInsert(nbnd, bndind, bndptr, higain);
 418
 419       edegrees = rinfo[higain].edegrees;
 420       edegrees[0] = edegrees[1] = 0;
 421       for (j=xadj[higain]; j<xadj[higain+1]; j++) {
 422         k = adjncy[j];
 423         if (where[k] == 2)
 424           rinfo[k].edegrees[to] -= vwgt[higain];
 425         else
 426           edegrees[where[k]] += vwgt[k];
 427       }
 428
 429       /* Push nodes out of the separator */
 430       for (j=mptr[nswaps]; j<mptr[nswaps+1]; j++) {
 431         k = mind[j];
 432         ASSERT(where[k] == 2);
 433         where[k] = from;
 434         INC_DEC(pwgts[from], pwgts[2], vwgt[k]);
 435         BNDDelete(nbnd, bndind, bndptr, k);
 436         for (jj=xadj[k]; jj<xadj[k+1]; jj++) {
 437           kk = adjncy[jj];
 438           if (where[kk] == 2)
 439             rinfo[kk].edegrees[from] += vwgt[k];
 440         }
 441       }
 442     }
 443
 444     ASSERT(mincut == pwgts[2]);
 445
 446     IFSET(ctrl->dbglvl, METIS_DBG_REFINE,
 447       printf("\tMinimum sep: %6"PRIDX" at %5"PRIDX", PWGTS: [%6"PRIDX" %6"PRIDX"], NBND: %6"PRIDX", QSIZE: %6"PRIDX"\n",
 448           mincut, mincutorder, pwgts[0], pwgts[1], nbnd, qsize));
 449
 450     graph->mincut = mincut;
 451     graph->nbnd   = nbnd;
 452
 453     if (pass%2 == 1 && (mincutorder == -1 || mincut >= initcut))
 454       break;
 455   }
 456
 457   rpqDestroy(queue);
 458
 459   WCOREPOP;
 460 }
 461
 462
 463 /*************************************************************************/
 464 /*! This function performs a node-based (two-sided) FM refinement that
 465     moves only nodes whose hmarker[] == -1. It is used by Parmetis. */
 466 /*************************************************************************/
 467 void FM_2WayNodeRefine2SidedP(ctrl_t *ctrl, graph_t *graph,
 468           idx_t *hmarker, real_t ubfactor, idx_t npasses)
 469 {
 470   idx_t i, ii, j, k, jj, kk, nvtxs, nbnd, nswaps, nmind;
 471   idx_t *xadj, *vwgt, *adjncy, *where, *pwgts, *edegrees, *bndind, *bndptr;
 472   idx_t *mptr, *mind, *moved, *swaps;
 473   rpq_t *queues[2];
 474   nrinfo_t *rinfo;
 475   idx_t higain, oldgain, mincut, initcut, mincutorder;
 476   idx_t pass, to, other, limit;
 477   idx_t badmaxpwgt, mindiff, newdiff;
 478   idx_t u[2], g[2];
 479
 480   WCOREPUSH;
 481
 482   nvtxs  = graph->nvtxs;
 483   xadj   = graph->xadj;
 484   adjncy = graph->adjncy;
 485   vwgt   = graph->vwgt;
 486
 487   bndind = graph->bndind;
 488   bndptr = graph->bndptr;
 489   where  = graph->where;
 490   pwgts  = graph->pwgts;
 491   rinfo  = graph->nrinfo;
 492
 493   queues[0] = rpqCreate(nvtxs);
 494   queues[1] = rpqCreate(nvtxs);
 495
 496   moved = iwspacemalloc(ctrl, nvtxs);
 497   swaps = iwspacemalloc(ctrl, nvtxs);
 498   mptr  = iwspacemalloc(ctrl, nvtxs+1);
 499   mind  = iwspacemalloc(ctrl, 2*nvtxs);
 500
 501   IFSET(ctrl->dbglvl, METIS_DBG_REFINE,
 502     printf("Partitions: [%6"PRIDX" %6"PRIDX"] Nv-Nb[%6"PRIDX" %6"PRIDX"]. ISep: %6"PRIDX"\n", pwgts[0], pwgts[1], graph->nvtxs, graph->nbnd, graph->mincut));
 503
 504   badmaxpwgt = (idx_t)(ubfactor*gk_max(pwgts[0], pwgts[1]));
 505
 506   for (pass=0; pass<npasses; pass++) {
 507     iset(nvtxs, -1, moved);
 508     rpqReset(queues[0]);
 509     rpqReset(queues[1]);
 510
 511     mincutorder = -1;
 512     initcut = mincut = graph->mincut;
 513     nbnd = graph->nbnd;
 514
 515     /* use the swaps array in place of the traditional perm array to save memory */
 516     irandArrayPermute(nbnd, swaps, nbnd, 1);
 517     for (ii=0; ii<nbnd; ii++) {
 518       i = bndind[swaps[ii]];
 519       ASSERT(where[i] == 2);
 520       if (hmarker[i] == -1) {
 521         rpqInsert(queues[0], i, vwgt[i]-rinfo[i].edegrees[1]);
 522         rpqInsert(queues[1], i, vwgt[i]-rinfo[i].edegrees[0]);
 523         moved[i] = -5;
 524       }
 525       else if (hmarker[i] != 2) {
 526         rpqInsert(queues[hmarker[i]], i, vwgt[i]-rinfo[i].edegrees[(hmarker[i]+1)%2]);
 527         moved[i] = -(10+hmarker[i]);
 528       }
 529     }
 530
 531     ASSERT(CheckNodeBnd(graph, nbnd));
 532     ASSERT(CheckNodePartitionParams(graph));
 533
 534     limit = nbnd;
 535
 536     /******************************************************
 537     * Get into the FM loop
 538     *******************************************************/
 539     mptr[0] = nmind = 0;
 540     mindiff = abs(pwgts[0]-pwgts[1]);
 541     to = (pwgts[0] < pwgts[1] ? 0 : 1);
 542     for (nswaps=0; nswaps<nvtxs; nswaps++) {
 543       u[0] = rpqSeeTopVal(queues[0]);
 544       u[1] = rpqSeeTopVal(queues[1]);
 545       if (u[0] != -1 && u[1] != -1) {
 546         g[0] = vwgt[u[0]]-rinfo[u[0]].edegrees[1];
 547         g[1] = vwgt[u[1]]-rinfo[u[1]].edegrees[0];
 548
 549         to = (g[0] > g[1] ? 0 : (g[0] < g[1] ? 1 : pass%2));
 550
 551         if (pwgts[to]+vwgt[u[to]] > badmaxpwgt)
 552           to = (to+1)%2;
 553       }
 554       else if (u[0] == -1 && u[1] == -1) {
 555         break;
 556       }
 557       else if (u[0] != -1 && pwgts[0]+vwgt[u[0]] <= badmaxpwgt) {
 558         to = 0;
 559       }
 560       else if (u[1] != -1 && pwgts[1]+vwgt[u[1]] <= badmaxpwgt) {
 561         to = 1;
 562       }
 563       else
 564         break;
 565
 566       other = (to+1)%2;
 567
 568       higain = rpqGetTop(queues[to]);
 569
 570       /* Delete its matching entry in the other queue */
 571       if (moved[higain] == -5)
 572         rpqDelete(queues[other], higain);
 573
 574       ASSERT(bndptr[higain] != -1);
 575
 576       /* The following check is to ensure we break out if there is a posibility
 577          of over-running the mind array.  */
 578       if (nmind + xadj[higain+1]-xadj[higain] >= 2*nvtxs-1)
 579         break;
 580
 581       pwgts[2] -= (vwgt[higain]-rinfo[higain].edegrees[other]);
 582
 583       newdiff = abs(pwgts[to]+vwgt[higain] - (pwgts[other]-rinfo[higain].edegrees[other]));
 584       if (pwgts[2] < mincut || (pwgts[2] == mincut && newdiff < mindiff)) {
 585         mincut      = pwgts[2];
 586         mincutorder = nswaps;
 587         mindiff     = newdiff;
 588       }
 589       else {
 590         if (nswaps - mincutorder > limit) {
 591           pwgts[2] += (vwgt[higain]-rinfo[higain].edegrees[other]);
 592           break; /* No further improvement, break out */
 593         }
 594       }
 595
 596       BNDDelete(nbnd, bndind, bndptr, higain);
 597       pwgts[to] += vwgt[higain];
 598       where[higain] = to;
 599       moved[higain] = nswaps;
 600       swaps[nswaps] = higain;
 601
 602
 603       /**********************************************************
 604       * Update the degrees of the affected nodes
 605       ***********************************************************/
 606       for (j=xadj[higain]; j<xadj[higain+1]; j++) {
 607         k = adjncy[j];
 608         if (where[k] == 2) { /* For the in-separator vertices modify their edegree[to] */
 609           oldgain = vwgt[k]-rinfo[k].edegrees[to];
 610           rinfo[k].edegrees[to] += vwgt[higain];
 611           if (moved[k] == -5 || moved[k] == -(10+other))
 612             rpqUpdate(queues[other], k, oldgain-vwgt[higain]);
 613         }
 614         else if (where[k] == other) { /* This vertex is pulled into the separator */
 615           ASSERTP(bndptr[k] == -1, ("%"PRIDX" %"PRIDX" %"PRIDX"\n", k, bndptr[k], where[k]));
 616           BNDInsert(nbnd, bndind, bndptr, k);
 617
 618           mind[nmind++] = k;  /* Keep track for rollback */
 619           where[k] = 2;
 620           pwgts[other] -= vwgt[k];
 621
 622           edegrees = rinfo[k].edegrees;
 623           edegrees[0] = edegrees[1] = 0;
 624           for (jj=xadj[k]; jj<xadj[k+1]; jj++) {
 625             kk = adjncy[jj];
 626             if (where[kk] != 2)
 627               edegrees[where[kk]] += vwgt[kk];
 628             else {
 629               oldgain = vwgt[kk]-rinfo[kk].edegrees[other];
 630               rinfo[kk].edegrees[other] -= vwgt[k];
 631               if (moved[kk] == -5 || moved[kk] == -(10+to))
 632                 rpqUpdate(queues[to], kk, oldgain+vwgt[k]);
 633             }
 634           }
 635
 636           /* Insert the new vertex into the priority queue (if it has not been moved). */
 637           if (moved[k] == -1 && (hmarker[k] == -1 || hmarker[k] == to)) {
 638             rpqInsert(queues[to], k, vwgt[k]-edegrees[other]);
 639             moved[k] = -(10+to);
 640           }
 641 #ifdef FULLMOVES  /* this does not work as well as the above partial one */
 642           if (moved[k] == -1) {
 643             if (hmarker[k] == -1) {
 644               rpqInsert(queues[0], k, vwgt[k]-edegrees[1]);
 645               rpqInsert(queues[1], k, vwgt[k]-edegrees[0]);
 646               moved[k] = -5;
 647             }
 648             else if (hmarker[k] != 2) {
 649               rpqInsert(queues[hmarker[k]], k, vwgt[k]-edegrees[(hmarker[k]+1)%2]);
 650               moved[k] = -(10+hmarker[k]);
 651             }
 652           }
 653 #endif
 654         }
 655       }
 656       mptr[nswaps+1] = nmind;
 657
 658       IFSET(ctrl->dbglvl, METIS_DBG_MOVEINFO,
 659             printf("Moved %6"PRIDX" to %3"PRIDX", Gain: %5"PRIDX" [%5"PRIDX"] "
 660                    "[%4"PRIDX" %4"PRIDX"] \t[%5"PRIDX" %5"PRIDX" %5"PRIDX"]\n",
 661                    higain, to, g[to], g[other], vwgt[u[to]], vwgt[u[other]],
 662                    pwgts[0], pwgts[1], pwgts[2]));
 663
 664     }
 665
 666
 667     /****************************************************************
 668     * Roll back computation
 669     *****************************************************************/
 670     for (nswaps--; nswaps>mincutorder; nswaps--) {
 671       higain = swaps[nswaps];
 672
 673       ASSERT(CheckNodePartitionParams(graph));
 674
 675       to = where[higain];
 676       other = (to+1)%2;
 677       INC_DEC(pwgts[2], pwgts[to], vwgt[higain]);
 678       where[higain] = 2;
 679       BNDInsert(nbnd, bndind, bndptr, higain);
 680
 681       edegrees = rinfo[higain].edegrees;
 682       edegrees[0] = edegrees[1] = 0;
 683       for (j=xadj[higain]; j<xadj[higain+1]; j++) {
 684         k = adjncy[j];
 685         if (where[k] == 2)
 686           rinfo[k].edegrees[to] -= vwgt[higain];
 687         else
 688           edegrees[where[k]] += vwgt[k];
 689       }
 690
 691       /* Push nodes out of the separator */
 692       for (j=mptr[nswaps]; j<mptr[nswaps+1]; j++) {
 693         k = mind[j];
 694         ASSERT(where[k] == 2);
 695         where[k] = other;
 696         INC_DEC(pwgts[other], pwgts[2], vwgt[k]);
 697         BNDDelete(nbnd, bndind, bndptr, k);
 698         for (jj=xadj[k]; jj<xadj[k+1]; jj++) {
 699           kk = adjncy[jj];
 700           if (where[kk] == 2)
 701             rinfo[kk].edegrees[other] += vwgt[k];
 702         }
 703       }
 704     }
 705
 706     ASSERT(mincut == pwgts[2]);
 707
 708     IFSET(ctrl->dbglvl, METIS_DBG_REFINE,
 709       printf("\tMinimum sep: %6"PRIDX" at %5"PRIDX", PWGTS: [%6"PRIDX" %6"PRIDX"], NBND: %6"PRIDX"\n", mincut, mincutorder, pwgts[0], pwgts[1], nbnd));
 710
 711     graph->mincut = mincut;
 712     graph->nbnd = nbnd;
 713
 714     if (mincutorder == -1 || mincut >= initcut)
 715       break;
 716   }
 717
 718   rpqDestroy(queues[0]);
 719   rpqDestroy(queues[1]);
 720
 721   WCOREPOP;
 722 }
 723