source: inundation/pymetis/metis-4.0/Lib/minitpart2.c @ 2051

Last change on this file since 2051 was 2051, checked in by jack, 19 years ago

Python interface to metis. Currently provides only the
METIS_PartMeshNodal function, since that is what is currently needed for partitioning.
Module name is metis.

File size: 10.7 KB
Line 
1/*
2 * Copyright 1997, Regents of the University of Minnesota
3 *
4 * minitpart2.c
5 *
6 * This file contains code that performs the initial partition of the
7 * coarsest graph
8 *
9 * Started 7/23/97
10 * George
11 *
12 * $Id: minitpart2.c,v 1.1 1998/11/27 17:59:23 karypis Exp $
13 *
14 */
15
16#include <metis.h>
17
18/*************************************************************************
19* This function computes the initial bisection of the coarsest graph
20**************************************************************************/
21void MocInit2WayPartition2(CtrlType *ctrl, GraphType *graph, float *tpwgts, float *ubvec) 
22{
23  int dbglvl;
24
25  dbglvl = ctrl->dbglvl;
26  IFSET(ctrl->dbglvl, DBG_REFINE, ctrl->dbglvl -= DBG_REFINE);
27  IFSET(ctrl->dbglvl, DBG_MOVEINFO, ctrl->dbglvl -= DBG_MOVEINFO);
28
29  IFSET(ctrl->dbglvl, DBG_TIME, starttimer(ctrl->InitPartTmr));
30
31  switch (ctrl->IType) {
32    case IPART_GGPKL:
33    case IPART_RANDOM:
34      MocGrowBisection2(ctrl, graph, tpwgts, ubvec);
35      break;
36    case 3:
37      MocGrowBisectionNew2(ctrl, graph, tpwgts, ubvec);
38      break;
39    default:
40      errexit("Unknown initial partition type: %d\n", ctrl->IType);
41  }
42
43  IFSET(ctrl->dbglvl, DBG_IPART, printf("Initial Cut: %d\n", graph->mincut));
44  IFSET(ctrl->dbglvl, DBG_TIME, stoptimer(ctrl->InitPartTmr));
45  ctrl->dbglvl = dbglvl;
46
47}
48
49
50
51
52/*************************************************************************
53* This function takes a graph and produces a bisection by using a region
54* growing algorithm. The resulting partition is returned in
55* graph->where
56**************************************************************************/
57void MocGrowBisection2(CtrlType *ctrl, GraphType *graph, float *tpwgts, float *ubvec)
58{
59  int i, j, k, nvtxs, ncon, from, bestcut, mincut, nbfs;
60  idxtype *bestwhere, *where;
61
62  nvtxs = graph->nvtxs;
63
64  MocAllocate2WayPartitionMemory(ctrl, graph);
65  where = graph->where;
66
67  bestwhere = idxmalloc(nvtxs, "BisectGraph: bestwhere");
68  nbfs = 2*(nvtxs <= ctrl->CoarsenTo ? SMALLNIPARTS : LARGENIPARTS);
69  bestcut = idxsum(graph->nedges, graph->adjwgt); 
70
71  for (; nbfs>0; nbfs--) {
72    idxset(nvtxs, 1, where);
73    where[RandomInRange(nvtxs)] = 0;
74
75    MocCompute2WayPartitionParams(ctrl, graph);
76
77    MocBalance2Way2(ctrl, graph, tpwgts, ubvec);
78
79    MocFM_2WayEdgeRefine2(ctrl, graph, tpwgts, ubvec, 4); 
80
81    MocBalance2Way2(ctrl, graph, tpwgts, ubvec);
82    MocFM_2WayEdgeRefine2(ctrl, graph, tpwgts, ubvec, 4); 
83
84    if (bestcut > graph->mincut) {
85      bestcut = graph->mincut;
86      idxcopy(nvtxs, where, bestwhere);
87      if (bestcut == 0)
88        break;
89    }
90  }
91
92  graph->mincut = bestcut;
93  idxcopy(nvtxs, bestwhere, where);
94
95  GKfree(&bestwhere, LTERM);
96}
97
98
99
100
101
102
103/*************************************************************************
104* This function takes a graph and produces a bisection by using a region
105* growing algorithm. The resulting partition is returned in
106* graph->where
107**************************************************************************/
108void MocGrowBisectionNew2(CtrlType *ctrl, GraphType *graph, float *tpwgts, float *ubvec)
109{
110  int i, j, k, nvtxs, ncon, from, bestcut, mincut, nbfs;
111  idxtype *bestwhere, *where;
112
113  nvtxs = graph->nvtxs;
114
115  MocAllocate2WayPartitionMemory(ctrl, graph);
116  where = graph->where;
117
118  bestwhere = idxmalloc(nvtxs, "BisectGraph: bestwhere");
119  nbfs = 2*(nvtxs <= ctrl->CoarsenTo ? SMALLNIPARTS : LARGENIPARTS);
120  bestcut = idxsum(graph->nedges, graph->adjwgt); 
121
122  for (; nbfs>0; nbfs--) {
123    idxset(nvtxs, 1, where);
124    where[RandomInRange(nvtxs)] = 0;
125
126    MocCompute2WayPartitionParams(ctrl, graph);
127
128    MocInit2WayBalance2(ctrl, graph, tpwgts, ubvec);
129
130    MocFM_2WayEdgeRefine2(ctrl, graph, tpwgts, ubvec, 4); 
131
132    if (bestcut > graph->mincut) {
133      bestcut = graph->mincut;
134      idxcopy(nvtxs, where, bestwhere);
135      if (bestcut == 0)
136        break;
137    }
138  }
139
140  graph->mincut = bestcut;
141  idxcopy(nvtxs, bestwhere, where);
142
143  GKfree(&bestwhere, LTERM);
144}
145
146
147
148/*************************************************************************
149* This function balances two partitions by moving the highest gain
150* (including negative gain) vertices to the other domain.
151* It is used only when tha unbalance is due to non contigous
152* subdomains. That is, the are no boundary vertices.
153* It moves vertices from the domain that is overweight to the one that
154* is underweight.
155**************************************************************************/
156void MocInit2WayBalance2(CtrlType *ctrl, GraphType *graph, float *tpwgts, float *ubvec)
157{
158  int i, ii, j, k, l, kwgt, nvtxs, nbnd, ncon, nswaps, from, to, pass, me, cnum, tmp, imin;
159  idxtype *xadj, *adjncy, *adjwgt, *where, *id, *ed, *bndptr, *bndind;
160  idxtype *moved, *perm, *qnum;
161  float *nvwgt, *npwgts, minwgt;
162  PQueueType parts[MAXNCON][2];
163  int higain, oldgain, mincut;
164
165  nvtxs = graph->nvtxs;
166  ncon = graph->ncon;
167  xadj = graph->xadj;
168  adjncy = graph->adjncy;
169  nvwgt = graph->nvwgt;
170  adjwgt = graph->adjwgt;
171  where = graph->where;
172  id = graph->id;
173  ed = graph->ed;
174  npwgts = graph->npwgts;
175  bndptr = graph->bndptr;
176  bndind = graph->bndind;
177
178  moved = idxwspacemalloc(ctrl, nvtxs);
179  perm = idxwspacemalloc(ctrl, nvtxs);
180  qnum = idxwspacemalloc(ctrl, nvtxs);
181
182  /* This is called for initial partitioning so we know from where to pick nodes */
183  from = 1;
184  to = (from+1)%2;
185
186  if (ctrl->dbglvl&DBG_REFINE) {
187    printf("Parts: [");
188    for (l=0; l<ncon; l++)
189      printf("(%.3f, %.3f) ", npwgts[l], npwgts[ncon+l]);
190    printf("] T[%.3f %.3f], Nv-Nb[%5d, %5d]. ICut: %6d, LB: %.3f [B]\n", tpwgts[0], tpwgts[1], graph->nvtxs, graph->nbnd, graph->mincut, ComputeLoadImbalance(ncon, 2, npwgts, tpwgts));
191  }
192
193  for (i=0; i<ncon; i++) {
194    PQueueInit(ctrl, &parts[i][0], nvtxs, PLUS_GAINSPAN+1);
195    PQueueInit(ctrl, &parts[i][1], nvtxs, PLUS_GAINSPAN+1);
196  }
197
198  idxset(nvtxs, -1, moved);
199
200  ASSERT(ComputeCut(graph, where) == graph->mincut);
201  ASSERT(CheckBnd(graph));
202  ASSERT(CheckGraph(graph));
203
204  /* Compute the queues in which each vertex will be assigned to */
205  for (i=0; i<nvtxs; i++)
206    qnum[i] = samax(ncon, nvwgt+i*ncon);
207
208  /* Insert the nodes of the proper partition in the appropriate priority queue */
209  RandomPermute(nvtxs, perm, 1);
210  for (ii=0; ii<nvtxs; ii++) {
211    i = perm[ii];
212    if (where[i] == from) {
213      if (ed[i] > 0)
214        PQueueInsert(&parts[qnum[i]][0], i, ed[i]-id[i]);
215      else
216        PQueueInsert(&parts[qnum[i]][1], i, ed[i]-id[i]);
217    }
218  }
219
220/*
221  for (i=0; i<ncon; i++)
222    printf("Queue #%d has %d %d\n", i, parts[i][0].nnodes, parts[i][1].nnodes);
223*/
224
225  /* Determine the termination criterion */
226  imin = 0;
227  for (i=1; i<ncon; i++) 
228    imin = (ubvec[i] < ubvec[imin] ? i : imin);
229  minwgt = .5/ubvec[imin];
230
231  mincut = graph->mincut;
232  nbnd = graph->nbnd;
233  for (nswaps=0; nswaps<nvtxs; nswaps++) {
234    /* Exit as soon as the minimum weight crossed over */
235    if (npwgts[to*ncon+imin] > minwgt) 
236      break;
237
238    if ((cnum = SelectQueueOneWay2(ncon, npwgts+to*ncon, parts, ubvec)) == -1)
239      break;
240
241    if ((higain = PQueueGetMax(&parts[cnum][0])) == -1)
242      higain = PQueueGetMax(&parts[cnum][1]);
243
244    mincut -= (ed[higain]-id[higain]);
245    saxpy(ncon, 1.0, nvwgt+higain*ncon, 1, npwgts+to*ncon, 1);
246    saxpy(ncon, -1.0, nvwgt+higain*ncon, 1, npwgts+from*ncon, 1);
247
248    where[higain] = to;
249    moved[higain] = nswaps;
250
251    if (ctrl->dbglvl&DBG_MOVEINFO) {
252      printf("Moved %6d from %d(%d). [%5d] %5d, NPwgts: ", higain, from, cnum, ed[higain]-id[higain], mincut);
253      for (l=0; l<ncon; l++) 
254        printf("(%.3f, %.3f) ", npwgts[l], npwgts[ncon+l]);
255      printf(", LB: %.3f\n", ComputeLoadImbalance(ncon, 2, npwgts, tpwgts));
256      if (ed[higain] == 0 && id[higain] > 0)
257        printf("\t Pulled from the interior!\n");
258    }
259
260
261    /**************************************************************
262    * Update the id[i]/ed[i] values of the affected nodes
263    ***************************************************************/
264    SWAP(id[higain], ed[higain], tmp);
265    if (ed[higain] == 0 && bndptr[higain] != -1 && xadj[higain] < xadj[higain+1]) 
266      BNDDelete(nbnd, bndind,  bndptr, higain);
267    if (ed[higain] > 0 && bndptr[higain] == -1)
268      BNDInsert(nbnd, bndind,  bndptr, higain);
269
270    for (j=xadj[higain]; j<xadj[higain+1]; j++) {
271      k = adjncy[j];
272      oldgain = ed[k]-id[k];
273
274      kwgt = (to == where[k] ? adjwgt[j] : -adjwgt[j]);
275      INC_DEC(id[k], ed[k], kwgt);
276
277      /* Update the queue position */
278      if (moved[k] == -1 && where[k] == from) {
279        if (ed[k] > 0 && bndptr[k] == -1) {  /* It moves in boundary */
280          PQueueDelete(&parts[qnum[k]][1], k, oldgain);
281          PQueueInsert(&parts[qnum[k]][0], k, ed[k]-id[k]);
282        }
283        else { /* It must be in the boundary already */
284          if (bndptr[k] == -1)
285            printf("What you thought was wrong!\n");
286          PQueueUpdate(&parts[qnum[k]][0], k, oldgain, ed[k]-id[k]);
287        }
288      }
289
290      /* Update its boundary information */
291      if (ed[k] == 0 && bndptr[k] != -1) 
292        BNDDelete(nbnd, bndind, bndptr, k);
293      else if (ed[k] > 0 && bndptr[k] == -1) 
294        BNDInsert(nbnd, bndind, bndptr, k);
295    }
296
297    ASSERTP(ComputeCut(graph, where) == mincut, ("%d != %d\n", ComputeCut(graph, where), mincut));
298
299  }
300
301  if (ctrl->dbglvl&DBG_REFINE) {
302    printf("\tMincut: %6d, NBND: %6d, NPwgts: ", mincut, nbnd);
303    for (l=0; l<ncon; l++)
304      printf("(%.3f, %.3f) ", npwgts[l], npwgts[ncon+l]);
305    printf(", LB: %.3f\n", ComputeLoadImbalance(ncon, 2, npwgts, tpwgts));
306  }
307
308  graph->mincut = mincut;
309  graph->nbnd = nbnd;
310
311  for (i=0; i<ncon; i++) {
312    PQueueFree(ctrl, &parts[i][0]);
313    PQueueFree(ctrl, &parts[i][1]);
314  }
315
316  ASSERT(ComputeCut(graph, where) == graph->mincut);
317  ASSERT(CheckBnd(graph));
318
319  idxwspacefree(ctrl, nvtxs);
320  idxwspacefree(ctrl, nvtxs);
321  idxwspacefree(ctrl, nvtxs);
322}
323
324
325
326/*************************************************************************
327* This function selects the partition number and the queue from which
328* we will move vertices out
329**************************************************************************/ 
330int SelectQueueOneWay2(int ncon, float *pto, PQueueType queues[MAXNCON][2], float *ubvec)
331{
332  int i, cnum=-1, imax, maxgain;
333  float max=0.0;
334  float twgt[MAXNCON];
335
336  for (i=0; i<ncon; i++) {
337    if (max < pto[i]) {
338      imax = i;
339      max = pto[i];
340    }
341  }
342  for (i=0; i<ncon; i++) 
343    twgt[i] = (max/(ubvec[imax]*ubvec[i]))/pto[i];
344  twgt[imax] = 0.0;
345
346  max = 0.0;
347  for (i=0; i<ncon; i++) {
348    if (max < twgt[i] && (PQueueGetSize(&queues[i][0]) > 0 || PQueueGetSize(&queues[i][1]) > 0)) {
349      max = twgt[i];
350      cnum = i;
351    }
352  }
353  if (max > 1)
354    return cnum;
355
356  /* optimize of cut */
357  maxgain = -10000000;
358  for (i=0; i<ncon; i++) {
359    if (PQueueGetSize(&queues[i][0]) > 0 && PQueueGetKey(&queues[i][0]) > maxgain) {
360      maxgain = PQueueGetKey(&queues[i][0]);
361      cnum = i;
362    }
363  }
364
365  return cnum;
366
367}
368
Note: See TracBrowser for help on using the repository browser.