Path: blob/devel/elmergrid/src/metis-5.1.0/libmetis/sfm.c
3206 views
/*1* Copyright 1997, Regents of the University of Minnesota2*3* sfm.c4*5* This file contains code that implementes an FM-based separator refinement6*7* Started 8/1/978* George9*10* $Id: sfm.c 10874 2011-10-17 23:13:00Z karypis $11*12*/1314#include "metislib.h"151617/*************************************************************************/18/*! This function performs a node-based FM refinement */19/**************************************************************************/20void FM_2WayNodeRefine2Sided(ctrl_t *ctrl, graph_t *graph, idx_t niter)21{22idx_t i, ii, j, k, jj, kk, nvtxs, nbnd, nswaps, nmind;23idx_t *xadj, *vwgt, *adjncy, *where, *pwgts, *edegrees, *bndind, *bndptr;24idx_t *mptr, *mind, *moved, *swaps;25rpq_t *queues[2];26nrinfo_t *rinfo;27idx_t higain, oldgain, mincut, initcut, mincutorder;28idx_t pass, to, other, limit;29idx_t badmaxpwgt, mindiff, newdiff;30idx_t u[2], g[2];31real_t mult;3233WCOREPUSH;3435nvtxs = graph->nvtxs;36xadj = graph->xadj;37adjncy = graph->adjncy;38vwgt = graph->vwgt;3940bndind = graph->bndind;41bndptr = graph->bndptr;42where = graph->where;43pwgts = graph->pwgts;44rinfo = graph->nrinfo;4546queues[0] = rpqCreate(nvtxs);47queues[1] = rpqCreate(nvtxs);4849moved = iwspacemalloc(ctrl, nvtxs);50swaps = iwspacemalloc(ctrl, nvtxs);51mptr = iwspacemalloc(ctrl, nvtxs+1);52mind = iwspacemalloc(ctrl, 2*nvtxs);5354mult = 0.5*ctrl->ubfactors[0];55badmaxpwgt = (idx_t)(mult*(pwgts[0]+pwgts[1]+pwgts[2]));5657IFSET(ctrl->dbglvl, METIS_DBG_REFINE,58printf("Partitions-N2: [%6"PRIDX" %6"PRIDX"] Nv-Nb[%6"PRIDX" %6"PRIDX"]. ISep: %6"PRIDX"\n", pwgts[0], pwgts[1], graph->nvtxs, graph->nbnd, graph->mincut));5960for (pass=0; pass<niter; pass++) {61iset(nvtxs, -1, moved);62rpqReset(queues[0]);63rpqReset(queues[1]);6465mincutorder = -1;66initcut = mincut = graph->mincut;67nbnd = graph->nbnd;6869/* use the swaps array in place of the traditional perm array to save memory */70irandArrayPermute(nbnd, swaps, nbnd, 1);71for (ii=0; ii<nbnd; ii++) {72i = bndind[swaps[ii]];73ASSERT(where[i] == 2);74rpqInsert(queues[0], i, vwgt[i]-rinfo[i].edegrees[1]);75rpqInsert(queues[1], i, vwgt[i]-rinfo[i].edegrees[0]);76}7778ASSERT(CheckNodeBnd(graph, nbnd));79ASSERT(CheckNodePartitionParams(graph));8081limit = (ctrl->compress ? gk_min(5*nbnd, 400) : gk_min(2*nbnd, 300));8283/******************************************************84* Get into the FM loop85*******************************************************/86mptr[0] = nmind = 0;87mindiff = iabs(pwgts[0]-pwgts[1]);88to = (pwgts[0] < pwgts[1] ? 0 : 1);89for (nswaps=0; nswaps<nvtxs; nswaps++) {90u[0] = rpqSeeTopVal(queues[0]);91u[1] = rpqSeeTopVal(queues[1]);92if (u[0] != -1 && u[1] != -1) {93g[0] = vwgt[u[0]]-rinfo[u[0]].edegrees[1];94g[1] = vwgt[u[1]]-rinfo[u[1]].edegrees[0];9596to = (g[0] > g[1] ? 0 : (g[0] < g[1] ? 1 : pass%2));9798if (pwgts[to]+vwgt[u[to]] > badmaxpwgt)99to = (to+1)%2;100}101else if (u[0] == -1 && u[1] == -1) {102break;103}104else if (u[0] != -1 && pwgts[0]+vwgt[u[0]] <= badmaxpwgt) {105to = 0;106}107else if (u[1] != -1 && pwgts[1]+vwgt[u[1]] <= badmaxpwgt) {108to = 1;109}110else111break;112113other = (to+1)%2;114115higain = rpqGetTop(queues[to]);116if (moved[higain] == -1) /* Delete if it was in the separator originally */117rpqDelete(queues[other], higain);118119ASSERT(bndptr[higain] != -1);120121/* The following check is to ensure we break out if there is a posibility122of over-running the mind array. */123if (nmind + xadj[higain+1]-xadj[higain] >= 2*nvtxs-1)124break;125126pwgts[2] -= (vwgt[higain]-rinfo[higain].edegrees[other]);127128newdiff = iabs(pwgts[to]+vwgt[higain] - (pwgts[other]-rinfo[higain].edegrees[other]));129if (pwgts[2] < mincut || (pwgts[2] == mincut && newdiff < mindiff)) {130mincut = pwgts[2];131mincutorder = nswaps;132mindiff = newdiff;133}134else {135if (nswaps - mincutorder > 2*limit ||136(nswaps - mincutorder > limit && pwgts[2] > 1.10*mincut)) {137pwgts[2] += (vwgt[higain]-rinfo[higain].edegrees[other]);138break; /* No further improvement, break out */139}140}141142BNDDelete(nbnd, bndind, bndptr, higain);143pwgts[to] += vwgt[higain];144where[higain] = to;145moved[higain] = nswaps;146swaps[nswaps] = higain;147148149/**********************************************************150* Update the degrees of the affected nodes151***********************************************************/152for (j=xadj[higain]; j<xadj[higain+1]; j++) {153k = adjncy[j];154if (where[k] == 2) { /* For the in-separator vertices modify their edegree[to] */155oldgain = vwgt[k]-rinfo[k].edegrees[to];156rinfo[k].edegrees[to] += vwgt[higain];157if (moved[k] == -1 || moved[k] == -(2+other))158rpqUpdate(queues[other], k, oldgain-vwgt[higain]);159}160else if (where[k] == other) { /* This vertex is pulled into the separator */161ASSERTP(bndptr[k] == -1, ("%"PRIDX" %"PRIDX" %"PRIDX"\n", k, bndptr[k], where[k]));162BNDInsert(nbnd, bndind, bndptr, k);163164mind[nmind++] = k; /* Keep track for rollback */165where[k] = 2;166pwgts[other] -= vwgt[k];167168edegrees = rinfo[k].edegrees;169edegrees[0] = edegrees[1] = 0;170for (jj=xadj[k]; jj<xadj[k+1]; jj++) {171kk = adjncy[jj];172if (where[kk] != 2)173edegrees[where[kk]] += vwgt[kk];174else {175oldgain = vwgt[kk]-rinfo[kk].edegrees[other];176rinfo[kk].edegrees[other] -= vwgt[k];177if (moved[kk] == -1 || moved[kk] == -(2+to))178rpqUpdate(queues[to], kk, oldgain+vwgt[k]);179}180}181182/* Insert the new vertex into the priority queue. Only one side! */183if (moved[k] == -1) {184rpqInsert(queues[to], k, vwgt[k]-edegrees[other]);185moved[k] = -(2+to);186}187}188}189mptr[nswaps+1] = nmind;190191IFSET(ctrl->dbglvl, METIS_DBG_MOVEINFO,192printf("Moved %6"PRIDX" to %3"PRIDX", Gain: %5"PRIDX" [%5"PRIDX"] [%4"PRIDX" %4"PRIDX"] \t[%5"PRIDX" %5"PRIDX" %5"PRIDX"]\n", higain, to, g[to], g[other], vwgt[u[to]], vwgt[u[other]], pwgts[0], pwgts[1], pwgts[2]));193194}195196197/****************************************************************198* Roll back computation199*****************************************************************/200for (nswaps--; nswaps>mincutorder; nswaps--) {201higain = swaps[nswaps];202203ASSERT(CheckNodePartitionParams(graph));204205to = where[higain];206other = (to+1)%2;207INC_DEC(pwgts[2], pwgts[to], vwgt[higain]);208where[higain] = 2;209BNDInsert(nbnd, bndind, bndptr, higain);210211edegrees = rinfo[higain].edegrees;212edegrees[0] = edegrees[1] = 0;213for (j=xadj[higain]; j<xadj[higain+1]; j++) {214k = adjncy[j];215if (where[k] == 2)216rinfo[k].edegrees[to] -= vwgt[higain];217else218edegrees[where[k]] += vwgt[k];219}220221/* Push nodes out of the separator */222for (j=mptr[nswaps]; j<mptr[nswaps+1]; j++) {223k = mind[j];224ASSERT(where[k] == 2);225where[k] = other;226INC_DEC(pwgts[other], pwgts[2], vwgt[k]);227BNDDelete(nbnd, bndind, bndptr, k);228for (jj=xadj[k]; jj<xadj[k+1]; jj++) {229kk = adjncy[jj];230if (where[kk] == 2)231rinfo[kk].edegrees[other] += vwgt[k];232}233}234}235236ASSERT(mincut == pwgts[2]);237238IFSET(ctrl->dbglvl, METIS_DBG_REFINE,239printf("\tMinimum sep: %6"PRIDX" at %5"PRIDX", PWGTS: [%6"PRIDX" %6"PRIDX"], NBND: %6"PRIDX"\n", mincut, mincutorder, pwgts[0], pwgts[1], nbnd));240241graph->mincut = mincut;242graph->nbnd = nbnd;243244if (mincutorder == -1 || mincut >= initcut)245break;246}247248rpqDestroy(queues[0]);249rpqDestroy(queues[1]);250251WCOREPOP;252}253254255/*************************************************************************/256/*! This function performs a node-based FM refinement.257Each refinement iteration is split into two sub-iterations.258In each sub-iteration only moves to one of the left/right partitions259is allowed; hence, it is one-sided.260*/261/**************************************************************************/262void FM_2WayNodeRefine1Sided(ctrl_t *ctrl, graph_t *graph, idx_t niter)263{264idx_t i, ii, j, k, jj, kk, nvtxs, nbnd, nswaps, nmind, iend;265idx_t *xadj, *vwgt, *adjncy, *where, *pwgts, *edegrees, *bndind, *bndptr;266idx_t *mptr, *mind, *swaps;267rpq_t *queue;268nrinfo_t *rinfo;269idx_t higain, mincut, initcut, mincutorder;270idx_t pass, to, other, limit;271idx_t badmaxpwgt, mindiff, newdiff;272real_t mult;273274WCOREPUSH;275276nvtxs = graph->nvtxs;277xadj = graph->xadj;278adjncy = graph->adjncy;279vwgt = graph->vwgt;280281bndind = graph->bndind;282bndptr = graph->bndptr;283where = graph->where;284pwgts = graph->pwgts;285rinfo = graph->nrinfo;286287queue = rpqCreate(nvtxs);288289swaps = iwspacemalloc(ctrl, nvtxs);290mptr = iwspacemalloc(ctrl, nvtxs+1);291mind = iwspacemalloc(ctrl, 2*nvtxs);292293mult = 0.5*ctrl->ubfactors[0];294badmaxpwgt = (idx_t)(mult*(pwgts[0]+pwgts[1]+pwgts[2]));295296IFSET(ctrl->dbglvl, METIS_DBG_REFINE,297printf("Partitions-N1: [%6"PRIDX" %6"PRIDX"] Nv-Nb[%6"PRIDX" %6"PRIDX"]. ISep: %6"PRIDX"\n", pwgts[0], pwgts[1], graph->nvtxs, graph->nbnd, graph->mincut));298299to = (pwgts[0] < pwgts[1] ? 1 : 0);300for (pass=0; pass<2*niter; pass++) { /* the 2*niter is for the two sides */301other = to;302to = (to+1)%2;303304rpqReset(queue);305306mincutorder = -1;307initcut = mincut = graph->mincut;308nbnd = graph->nbnd;309310/* use the swaps array in place of the traditional perm array to save memory */311irandArrayPermute(nbnd, swaps, nbnd, 1);312for (ii=0; ii<nbnd; ii++) {313i = bndind[swaps[ii]];314ASSERT(where[i] == 2);315rpqInsert(queue, i, vwgt[i]-rinfo[i].edegrees[other]);316}317318ASSERT(CheckNodeBnd(graph, nbnd));319ASSERT(CheckNodePartitionParams(graph));320321limit = (ctrl->compress ? gk_min(5*nbnd, 500) : gk_min(3*nbnd, 300));322323/******************************************************324* Get into the FM loop325*******************************************************/326IFSET(ctrl->dbglvl, METIS_DBG_TIME, gk_startcputimer(ctrl->Aux3Tmr));327mptr[0] = nmind = 0;328mindiff = iabs(pwgts[0]-pwgts[1]);329for (nswaps=0; nswaps<nvtxs; nswaps++) {330if ((higain = rpqGetTop(queue)) == -1)331break;332333ASSERT(bndptr[higain] != -1);334335/* The following check is to ensure we break out if there is a posibility336of over-running the mind array. */337if (nmind + xadj[higain+1]-xadj[higain] >= 2*nvtxs-1)338break;339340if (pwgts[to]+vwgt[higain] > badmaxpwgt)341break; /* No point going any further. Balance will be bad */342343pwgts[2] -= (vwgt[higain]-rinfo[higain].edegrees[other]);344345newdiff = iabs(pwgts[to]+vwgt[higain] - (pwgts[other]-rinfo[higain].edegrees[other]));346if (pwgts[2] < mincut || (pwgts[2] == mincut && newdiff < mindiff)) {347mincut = pwgts[2];348mincutorder = nswaps;349mindiff = newdiff;350}351else {352if (nswaps - mincutorder > 3*limit ||353(nswaps - mincutorder > limit && pwgts[2] > 1.10*mincut)) {354pwgts[2] += (vwgt[higain]-rinfo[higain].edegrees[other]);355break; /* No further improvement, break out */356}357}358359BNDDelete(nbnd, bndind, bndptr, higain);360pwgts[to] += vwgt[higain];361where[higain] = to;362swaps[nswaps] = higain;363364365/**********************************************************366* Update the degrees of the affected nodes367***********************************************************/368IFSET(ctrl->dbglvl, METIS_DBG_TIME, gk_startcputimer(ctrl->Aux1Tmr));369for (j=xadj[higain]; j<xadj[higain+1]; j++) {370k = adjncy[j];371372if (where[k] == 2) { /* For the in-separator vertices modify their edegree[to] */373rinfo[k].edegrees[to] += vwgt[higain];374}375else if (where[k] == other) { /* This vertex is pulled into the separator */376ASSERTP(bndptr[k] == -1, ("%"PRIDX" %"PRIDX" %"PRIDX"\n", k, bndptr[k], where[k]));377BNDInsert(nbnd, bndind, bndptr, k);378379mind[nmind++] = k; /* Keep track for rollback */380where[k] = 2;381pwgts[other] -= vwgt[k];382383edegrees = rinfo[k].edegrees;384edegrees[0] = edegrees[1] = 0;385for (jj=xadj[k], iend=xadj[k+1]; jj<iend; jj++) {386kk = adjncy[jj];387if (where[kk] != 2)388edegrees[where[kk]] += vwgt[kk];389else {390rinfo[kk].edegrees[other] -= vwgt[k];391392/* Since the moves are one-sided this vertex has not been moved yet */393rpqUpdate(queue, kk, vwgt[kk]-rinfo[kk].edegrees[other]);394}395}396397/* Insert the new vertex into the priority queue. Safe due to one-sided moves */398rpqInsert(queue, k, vwgt[k]-edegrees[other]);399}400}401mptr[nswaps+1] = nmind;402IFSET(ctrl->dbglvl, METIS_DBG_TIME, gk_stopcputimer(ctrl->Aux1Tmr));403404405IFSET(ctrl->dbglvl, METIS_DBG_MOVEINFO,406printf("Moved %6"PRIDX" to %3"PRIDX", Gain: %5"PRIDX" [%5"PRIDX"] \t[%5"PRIDX" %5"PRIDX" %5"PRIDX"] [%3"PRIDX" %2"PRIDX"]\n",407higain, to, (vwgt[higain]-rinfo[higain].edegrees[other]), vwgt[higain],408pwgts[0], pwgts[1], pwgts[2], nswaps, limit));409}410IFSET(ctrl->dbglvl, METIS_DBG_TIME, gk_stopcputimer(ctrl->Aux3Tmr));411412413/****************************************************************414* Roll back computation415*****************************************************************/416IFSET(ctrl->dbglvl, METIS_DBG_TIME, gk_startcputimer(ctrl->Aux2Tmr));417for (nswaps--; nswaps>mincutorder; nswaps--) {418higain = swaps[nswaps];419420ASSERT(CheckNodePartitionParams(graph));421ASSERT(where[higain] == to);422423INC_DEC(pwgts[2], pwgts[to], vwgt[higain]);424where[higain] = 2;425BNDInsert(nbnd, bndind, bndptr, higain);426427edegrees = rinfo[higain].edegrees;428edegrees[0] = edegrees[1] = 0;429for (j=xadj[higain]; j<xadj[higain+1]; j++) {430k = adjncy[j];431if (where[k] == 2)432rinfo[k].edegrees[to] -= vwgt[higain];433else434edegrees[where[k]] += vwgt[k];435}436437/* Push nodes out of the separator */438for (j=mptr[nswaps]; j<mptr[nswaps+1]; j++) {439k = mind[j];440ASSERT(where[k] == 2);441where[k] = other;442INC_DEC(pwgts[other], pwgts[2], vwgt[k]);443BNDDelete(nbnd, bndind, bndptr, k);444for (jj=xadj[k], iend=xadj[k+1]; jj<iend; jj++) {445kk = adjncy[jj];446if (where[kk] == 2)447rinfo[kk].edegrees[other] += vwgt[k];448}449}450}451IFSET(ctrl->dbglvl, METIS_DBG_TIME, gk_stopcputimer(ctrl->Aux2Tmr));452453ASSERT(mincut == pwgts[2]);454455IFSET(ctrl->dbglvl, METIS_DBG_REFINE,456printf("\tMinimum sep: %6"PRIDX" at %5"PRIDX", PWGTS: [%6"PRIDX" %6"PRIDX"], NBND: %6"PRIDX"\n", mincut, mincutorder, pwgts[0], pwgts[1], nbnd));457458graph->mincut = mincut;459graph->nbnd = nbnd;460461if (pass%2 == 1 && (mincutorder == -1 || mincut >= initcut))462break;463}464465rpqDestroy(queue);466467WCOREPOP;468}469470471/*************************************************************************/472/*! This function balances the left/right partitions of a separator473tri-section */474/*************************************************************************/475void FM_2WayNodeBalance(ctrl_t *ctrl, graph_t *graph)476{477idx_t i, ii, j, k, jj, kk, nvtxs, nbnd, nswaps, gain;478idx_t badmaxpwgt, higain, oldgain, pass, to, other;479idx_t *xadj, *vwgt, *adjncy, *where, *pwgts, *edegrees, *bndind, *bndptr;480idx_t *perm, *moved;481rpq_t *queue;482nrinfo_t *rinfo;483real_t mult;484485nvtxs = graph->nvtxs;486xadj = graph->xadj;487adjncy = graph->adjncy;488vwgt = graph->vwgt;489490bndind = graph->bndind;491bndptr = graph->bndptr;492where = graph->where;493pwgts = graph->pwgts;494rinfo = graph->nrinfo;495496mult = 0.5*ctrl->ubfactors[0];497498badmaxpwgt = (idx_t)(mult*(pwgts[0]+pwgts[1]));499if (gk_max(pwgts[0], pwgts[1]) < badmaxpwgt)500return;501if (iabs(pwgts[0]-pwgts[1]) < 3*graph->tvwgt[0]/nvtxs)502return;503504WCOREPUSH;505506to = (pwgts[0] < pwgts[1] ? 0 : 1);507other = (to+1)%2;508509queue = rpqCreate(nvtxs);510511perm = iwspacemalloc(ctrl, nvtxs);512moved = iset(nvtxs, -1, iwspacemalloc(ctrl, nvtxs));513514IFSET(ctrl->dbglvl, METIS_DBG_REFINE,515printf("Partitions: [%6"PRIDX" %6"PRIDX"] Nv-Nb[%6"PRIDX" %6"PRIDX"]. ISep: %6"PRIDX" [B]\n", pwgts[0], pwgts[1], graph->nvtxs, graph->nbnd, graph->mincut));516517nbnd = graph->nbnd;518irandArrayPermute(nbnd, perm, nbnd, 1);519for (ii=0; ii<nbnd; ii++) {520i = bndind[perm[ii]];521ASSERT(where[i] == 2);522rpqInsert(queue, i, vwgt[i]-rinfo[i].edegrees[other]);523}524525ASSERT(CheckNodeBnd(graph, nbnd));526ASSERT(CheckNodePartitionParams(graph));527528/******************************************************529* Get into the FM loop530*******************************************************/531for (nswaps=0; nswaps<nvtxs; nswaps++) {532if ((higain = rpqGetTop(queue)) == -1)533break;534535moved[higain] = 1;536537gain = vwgt[higain]-rinfo[higain].edegrees[other];538badmaxpwgt = (idx_t)(mult*(pwgts[0]+pwgts[1]));539540/* break if other is now underwight */541if (pwgts[to] > pwgts[other])542break;543544/* break if balance is achieved and no +ve or zero gain */545if (gain < 0 && pwgts[other] < badmaxpwgt)546break;547548/* skip this vertex if it will violate balance on the other side */549if (pwgts[to]+vwgt[higain] > badmaxpwgt)550continue;551552ASSERT(bndptr[higain] != -1);553554pwgts[2] -= gain;555556BNDDelete(nbnd, bndind, bndptr, higain);557pwgts[to] += vwgt[higain];558where[higain] = to;559560IFSET(ctrl->dbglvl, METIS_DBG_MOVEINFO,561printf("Moved %6"PRIDX" to %3"PRIDX", Gain: %3"PRIDX", \t[%5"PRIDX" %5"PRIDX" %5"PRIDX"]\n", higain, to, vwgt[higain]-rinfo[higain].edegrees[other], pwgts[0], pwgts[1], pwgts[2]));562563564/**********************************************************565* Update the degrees of the affected nodes566***********************************************************/567for (j=xadj[higain]; j<xadj[higain+1]; j++) {568k = adjncy[j];569if (where[k] == 2) { /* For the in-separator vertices modify their edegree[to] */570rinfo[k].edegrees[to] += vwgt[higain];571}572else if (where[k] == other) { /* This vertex is pulled into the separator */573ASSERTP(bndptr[k] == -1, ("%"PRIDX" %"PRIDX" %"PRIDX"\n", k, bndptr[k], where[k]));574BNDInsert(nbnd, bndind, bndptr, k);575576where[k] = 2;577pwgts[other] -= vwgt[k];578579edegrees = rinfo[k].edegrees;580edegrees[0] = edegrees[1] = 0;581for (jj=xadj[k]; jj<xadj[k+1]; jj++) {582kk = adjncy[jj];583if (where[kk] != 2)584edegrees[where[kk]] += vwgt[kk];585else {586ASSERT(bndptr[kk] != -1);587oldgain = vwgt[kk]-rinfo[kk].edegrees[other];588rinfo[kk].edegrees[other] -= vwgt[k];589590if (moved[kk] == -1)591rpqUpdate(queue, kk, oldgain+vwgt[k]);592}593}594595/* Insert the new vertex into the priority queue */596rpqInsert(queue, k, vwgt[k]-edegrees[other]);597}598}599}600601IFSET(ctrl->dbglvl, METIS_DBG_REFINE,602printf("\tBalanced sep: %6"PRIDX" at %4"PRIDX", PWGTS: [%6"PRIDX" %6"PRIDX"], NBND: %6"PRIDX"\n", pwgts[2], nswaps, pwgts[0], pwgts[1], nbnd));603604graph->mincut = pwgts[2];605graph->nbnd = nbnd;606607rpqDestroy(queue);608609WCOREPOP;610}611612613614