1 /**
2 \file
3 \brief This file contains various routines for dealing with options and ctrl_t.
4
5 \date Started 5/12/2011
6 \author George
7 \author Copyright 1997-2011, Regents of the University of Minnesota
8 */
9
10 #include "metislib.h"
11
12
13 /*************************************************************************/
14 /*! This function creates and sets the run parameters (ctrl_t) */
15 /*************************************************************************/
SetupCtrl(moptype_et optype,idx_t * options,idx_t ncon,idx_t nparts,real_t * tpwgts,real_t * ubvec)16 ctrl_t *SetupCtrl(moptype_et optype, idx_t *options, idx_t ncon, idx_t nparts,
17 real_t *tpwgts, real_t *ubvec)
18 {
19 idx_t i, j;
20 ctrl_t *ctrl;
21
22 ctrl = (ctrl_t *)gk_malloc(sizeof(ctrl_t), "SetupCtrl: ctrl");
23
24 memset((void *)ctrl, 0, sizeof(ctrl_t));
25
26 switch (optype) {
27 case METIS_OP_PMETIS:
28 ctrl->objtype = GETOPTION(options, METIS_OPTION_OBJTYPE, METIS_OBJTYPE_CUT);
29 ctrl->rtype = METIS_RTYPE_FM;
30 ctrl->ncuts = GETOPTION(options, METIS_OPTION_NCUTS, 1);
31 ctrl->niter = GETOPTION(options, METIS_OPTION_NITER, 10);
32
33 if (ncon == 1) {
34 ctrl->iptype = GETOPTION(options, METIS_OPTION_IPTYPE, METIS_IPTYPE_GROW);
35 ctrl->ufactor = GETOPTION(options, METIS_OPTION_UFACTOR, PMETIS_DEFAULT_UFACTOR);
36 ctrl->CoarsenTo = 20;
37 }
38 else {
39 ctrl->iptype = GETOPTION(options, METIS_OPTION_IPTYPE, METIS_IPTYPE_RANDOM);
40 ctrl->ufactor = GETOPTION(options, METIS_OPTION_UFACTOR, MCPMETIS_DEFAULT_UFACTOR);
41 ctrl->CoarsenTo = 100;
42 }
43
44 break;
45
46
47 case METIS_OP_KMETIS:
48 ctrl->objtype = GETOPTION(options, METIS_OPTION_OBJTYPE, METIS_OBJTYPE_CUT);
49 ctrl->iptype = METIS_IPTYPE_METISRB;
50 ctrl->rtype = METIS_RTYPE_GREEDY;
51 ctrl->ncuts = GETOPTION(options, METIS_OPTION_NCUTS, 1);
52 ctrl->niter = GETOPTION(options, METIS_OPTION_NITER, 10);
53 ctrl->ufactor = GETOPTION(options, METIS_OPTION_UFACTOR, KMETIS_DEFAULT_UFACTOR);
54 ctrl->minconn = GETOPTION(options, METIS_OPTION_MINCONN, 0);
55 ctrl->contig = GETOPTION(options, METIS_OPTION_CONTIG, 0);
56 break;
57
58
59 case METIS_OP_OMETIS:
60 ctrl->objtype = GETOPTION(options, METIS_OPTION_OBJTYPE, METIS_OBJTYPE_NODE);
61 ctrl->rtype = GETOPTION(options, METIS_OPTION_RTYPE, METIS_RTYPE_SEP1SIDED);
62 ctrl->iptype = GETOPTION(options, METIS_OPTION_IPTYPE, METIS_IPTYPE_EDGE);
63 ctrl->nseps = GETOPTION(options, METIS_OPTION_NSEPS, 1);
64 ctrl->niter = GETOPTION(options, METIS_OPTION_NITER, 10);
65 ctrl->ufactor = GETOPTION(options, METIS_OPTION_UFACTOR, OMETIS_DEFAULT_UFACTOR);
66 ctrl->compress = GETOPTION(options, METIS_OPTION_COMPRESS, 1);
67 ctrl->ccorder = GETOPTION(options, METIS_OPTION_CCORDER, 0);
68 ctrl->pfactor = 0.1*GETOPTION(options, METIS_OPTION_PFACTOR, 0);
69
70 ctrl->CoarsenTo = 100;
71 break;
72
73 default:
74 gk_errexit(SIGERR, "Unknown optype of %d\n", optype);
75 }
76
77 /* common options */
78 ctrl->ctype = GETOPTION(options, METIS_OPTION_CTYPE, METIS_CTYPE_SHEM);
79 ctrl->no2hop = GETOPTION(options, METIS_OPTION_NO2HOP, 0);
80 ctrl->seed = GETOPTION(options, METIS_OPTION_SEED, -1);
81 ctrl->dbglvl = GETOPTION(options, METIS_OPTION_DBGLVL, 0);
82 ctrl->numflag = GETOPTION(options, METIS_OPTION_NUMBERING, 0);
83
84 /* set non-option information */
85 ctrl->optype = optype;
86 ctrl->ncon = ncon;
87 ctrl->nparts = nparts;
88 ctrl->maxvwgt = ismalloc(ncon, 0, "SetupCtrl: maxvwgt");
89
90 /* setup the target partition weights */
91 if (ctrl->optype != METIS_OP_OMETIS) {
92 ctrl->tpwgts = rmalloc(nparts*ncon, "SetupCtrl: ctrl->tpwgts");
93 if (tpwgts) {
94 rcopy(nparts*ncon, tpwgts, ctrl->tpwgts);
95 }
96 else {
97 for (i=0; i<nparts; i++) {
98 for (j=0; j<ncon; j++)
99 ctrl->tpwgts[i*ncon+j] = 1.0/nparts;
100 }
101 }
102 }
103 else { /* METIS_OP_OMETIS */
104 /* this is required to allow the pijbm to be defined properly for
105 the edge-based refinement during initial partitioning */
106 ctrl->tpwgts = rsmalloc(2, .5, "SetupCtrl: ctrl->tpwgts");
107 }
108
109
110 /* setup the ubfactors */
111 ctrl->ubfactors = rsmalloc(ctrl->ncon, I2RUBFACTOR(ctrl->ufactor), "SetupCtrl: ubfactors");
112 if (ubvec)
113 rcopy(ctrl->ncon, ubvec, ctrl->ubfactors);
114 for (i=0; i<ctrl->ncon; i++)
115 ctrl->ubfactors[i] += 0.0000499;
116
117 /* Allocate memory for balance multipliers.
118 Note that for PMETIS/OMETIS routines the memory allocated is more
119 than required as balance multipliers for 2 parts is sufficient. */
120 ctrl->pijbm = rmalloc(nparts*ncon, "SetupCtrl: ctrl->pijbm");
121
122 InitRandom(ctrl->seed);
123
124 IFSET(ctrl->dbglvl, METIS_DBG_INFO, PrintCtrl(ctrl));
125
126 if (!CheckParams(ctrl)) {
127 FreeCtrl(&ctrl);
128 return NULL;
129 }
130 else {
131 return ctrl;
132 }
133 }
134
135
136 /*************************************************************************/
137 /*! Computes the per-partition/constraint balance multipliers */
138 /*************************************************************************/
SetupKWayBalMultipliers(ctrl_t * ctrl,graph_t * graph)139 void SetupKWayBalMultipliers(ctrl_t *ctrl, graph_t *graph)
140 {
141 idx_t i, j;
142
143 for (i=0; i<ctrl->nparts; i++) {
144 for (j=0; j<graph->ncon; j++)
145 ctrl->pijbm[i*graph->ncon+j] = graph->invtvwgt[j]/ctrl->tpwgts[i*graph->ncon+j];
146 }
147 }
148
149
150 /*************************************************************************/
151 /*! Computes the per-partition/constraint balance multipliers */
152 /*************************************************************************/
Setup2WayBalMultipliers(ctrl_t * ctrl,graph_t * graph,real_t * tpwgts)153 void Setup2WayBalMultipliers(ctrl_t *ctrl, graph_t *graph, real_t *tpwgts)
154 {
155 idx_t i, j;
156
157 for (i=0; i<2; i++) {
158 for (j=0; j<graph->ncon; j++)
159 ctrl->pijbm[i*graph->ncon+j] = graph->invtvwgt[j]/tpwgts[i*graph->ncon+j];
160 }
161 }
162
163
164 /*************************************************************************/
165 /*! This function prints the various control fields */
166 /*************************************************************************/
PrintCtrl(ctrl_t * ctrl)167 void PrintCtrl(ctrl_t *ctrl)
168 {
169 idx_t i, j, modnum;
170
171 printf(" Runtime parameters:\n");
172
173 printf(" Objective type: ");
174 switch (ctrl->objtype) {
175 case METIS_OBJTYPE_CUT:
176 printf("METIS_OBJTYPE_CUT\n");
177 break;
178 case METIS_OBJTYPE_VOL:
179 printf("METIS_OBJTYPE_VOL\n");
180 break;
181 case METIS_OBJTYPE_NODE:
182 printf("METIS_OBJTYPE_NODE\n");
183 break;
184 default:
185 printf("Unknown!\n");
186 }
187
188 printf(" Coarsening type: ");
189 switch (ctrl->ctype) {
190 case METIS_CTYPE_RM:
191 printf("METIS_CTYPE_RM\n");
192 break;
193 case METIS_CTYPE_SHEM:
194 printf("METIS_CTYPE_SHEM\n");
195 break;
196 default:
197 printf("Unknown!\n");
198 }
199
200 printf(" Initial partitioning type: ");
201 switch (ctrl->iptype) {
202 case METIS_IPTYPE_GROW:
203 printf("METIS_IPTYPE_GROW\n");
204 break;
205 case METIS_IPTYPE_RANDOM:
206 printf("METIS_IPTYPE_RANDOM\n");
207 break;
208 case METIS_IPTYPE_EDGE:
209 printf("METIS_IPTYPE_EDGE\n");
210 break;
211 case METIS_IPTYPE_NODE:
212 printf("METIS_IPTYPE_NODE\n");
213 break;
214 case METIS_IPTYPE_METISRB:
215 printf("METIS_IPTYPE_METISRB\n");
216 break;
217 default:
218 printf("Unknown!\n");
219 }
220
221 printf(" Refinement type: ");
222 switch (ctrl->rtype) {
223 case METIS_RTYPE_FM:
224 printf("METIS_RTYPE_FM\n");
225 break;
226 case METIS_RTYPE_GREEDY:
227 printf("METIS_RTYPE_GREEDY\n");
228 break;
229 case METIS_RTYPE_SEP2SIDED:
230 printf("METIS_RTYPE_SEP2SIDED\n");
231 break;
232 case METIS_RTYPE_SEP1SIDED:
233 printf("METIS_RTYPE_SEP1SIDED\n");
234 break;
235 default:
236 printf("Unknown!\n");
237 }
238
239 printf(" Perform a 2-hop matching: %s\n", (ctrl->no2hop ? "Yes" : "No"));
240
241 printf(" Number of balancing constraints: %"PRIDX"\n", ctrl->ncon);
242 printf(" Number of refinement iterations: %"PRIDX"\n", ctrl->niter);
243 printf(" Random number seed: %"PRIDX"\n", ctrl->seed);
244
245 if (ctrl->optype == METIS_OP_OMETIS) {
246 printf(" Number of separators: %"PRIDX"\n", ctrl->nseps);
247 printf(" Compress graph prior to ordering: %s\n", (ctrl->compress ? "Yes" : "No"));
248 printf(" Detect & order connected components separately: %s\n", (ctrl->ccorder ? "Yes" : "No"));
249 printf(" Prunning factor for high degree vertices: %"PRREAL"\n", ctrl->pfactor);
250 }
251 else {
252 printf(" Number of partitions: %"PRIDX"\n", ctrl->nparts);
253 printf(" Number of cuts: %"PRIDX"\n", ctrl->ncuts);
254 printf(" User-supplied ufactor: %"PRIDX"\n", ctrl->ufactor);
255
256 if (ctrl->optype == METIS_OP_KMETIS) {
257 printf(" Minimize connectivity: %s\n", (ctrl->minconn ? "Yes" : "No"));
258 printf(" Create contigous partitions: %s\n", (ctrl->contig ? "Yes" : "No"));
259 }
260
261 modnum = (ctrl->ncon==1 ? 5 : (ctrl->ncon==2 ? 3 : (ctrl->ncon==3 ? 2 : 1)));
262 printf(" Target partition weights: ");
263 for (i=0; i<ctrl->nparts; i++) {
264 if (i%modnum == 0)
265 printf("\n ");
266 printf("%4"PRIDX"=[", i);
267 for (j=0; j<ctrl->ncon; j++)
268 printf("%s%.2e", (j==0 ? "" : " "), (double)ctrl->tpwgts[i*ctrl->ncon+j]);
269 printf("]");
270 }
271 printf("\n");
272 }
273
274 printf(" Allowed maximum load imbalance: ");
275 for (i=0; i<ctrl->ncon; i++)
276 printf("%.3"PRREAL" ", ctrl->ubfactors[i]);
277 printf("\n");
278
279 printf("\n");
280 }
281
282
283 /*************************************************************************/
284 /*! This function checks the validity of user-supplied parameters */
285 /*************************************************************************/
CheckParams(ctrl_t * ctrl)286 int CheckParams(ctrl_t *ctrl)
287 {
288 idx_t i, j;
289 real_t sum;
290 mdbglvl_et dbglvl=METIS_DBG_INFO;
291
292 switch (ctrl->optype) {
293 case METIS_OP_PMETIS:
294 if (ctrl->objtype != METIS_OBJTYPE_CUT) {
295 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect objective type.\n"));
296 return 0;
297 }
298 if (ctrl->ctype != METIS_CTYPE_RM && ctrl->ctype != METIS_CTYPE_SHEM) {
299 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect coarsening scheme.\n"));
300 return 0;
301 }
302 if (ctrl->iptype != METIS_IPTYPE_GROW && ctrl->iptype != METIS_IPTYPE_RANDOM) {
303 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect initial partitioning scheme.\n"));
304 return 0;
305 }
306 if (ctrl->rtype != METIS_RTYPE_FM) {
307 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect refinement scheme.\n"));
308 return 0;
309 }
310 if (ctrl->ncuts <= 0) {
311 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect ncuts.\n"));
312 return 0;
313 }
314 if (ctrl->niter <= 0) {
315 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect niter.\n"));
316 return 0;
317 }
318 if (ctrl->ufactor <= 0) {
319 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect ufactor.\n"));
320 return 0;
321 }
322 if (ctrl->numflag != 0 && ctrl->numflag != 1) {
323 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect numflag.\n"));
324 return 0;
325 }
326 if (ctrl->nparts <= 0) {
327 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect nparts.\n"));
328 return 0;
329 }
330 if (ctrl->ncon <= 0) {
331 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect ncon.\n"));
332 return 0;
333 }
334
335 for (i=0; i<ctrl->ncon; i++) {
336 sum = rsum(ctrl->nparts, ctrl->tpwgts+i, ctrl->ncon);
337 if (sum < 0.99 || sum > 1.01) {
338 IFSET(dbglvl, METIS_DBG_INFO,
339 printf("Input Error: Incorrect sum of %"PRREAL" for tpwgts for constraint %"PRIDX".\n", sum, i));
340 return 0;
341 }
342 }
343 for (i=0; i<ctrl->ncon; i++) {
344 for (j=0; j<ctrl->nparts; j++) {
345 if (ctrl->tpwgts[j*ctrl->ncon+i] <= 0.0) {
346 IFSET(dbglvl, METIS_DBG_INFO,
347 printf("Input Error: Incorrect tpwgts for partition %"PRIDX" and constraint %"PRIDX".\n", j, i));
348 return 0;
349 }
350 }
351 }
352
353 for (i=0; i<ctrl->ncon; i++) {
354 if (ctrl->ubfactors[i] <= 1.0) {
355 IFSET(dbglvl, METIS_DBG_INFO,
356 printf("Input Error: Incorrect ubfactor for constraint %"PRIDX".\n", i));
357 return 0;
358 }
359 }
360
361 break;
362
363 case METIS_OP_KMETIS:
364 if (ctrl->objtype != METIS_OBJTYPE_CUT && ctrl->objtype != METIS_OBJTYPE_VOL) {
365 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect objective type.\n"));
366 return 0;
367 }
368 if (ctrl->ctype != METIS_CTYPE_RM && ctrl->ctype != METIS_CTYPE_SHEM) {
369 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect coarsening scheme.\n"));
370 return 0;
371 }
372 if (ctrl->iptype != METIS_IPTYPE_METISRB) {
373 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect initial partitioning scheme.\n"));
374 return 0;
375 }
376 if (ctrl->rtype != METIS_RTYPE_GREEDY) {
377 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect refinement scheme.\n"));
378 return 0;
379 }
380 if (ctrl->ncuts <= 0) {
381 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect ncuts.\n"));
382 return 0;
383 }
384 if (ctrl->niter <= 0) {
385 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect niter.\n"));
386 return 0;
387 }
388 if (ctrl->ufactor <= 0) {
389 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect ufactor.\n"));
390 return 0;
391 }
392 if (ctrl->numflag != 0 && ctrl->numflag != 1) {
393 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect numflag.\n"));
394 return 0;
395 }
396 if (ctrl->nparts <= 0) {
397 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect nparts.\n"));
398 return 0;
399 }
400 if (ctrl->ncon <= 0) {
401 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect ncon.\n"));
402 return 0;
403 }
404 if (ctrl->contig != 0 && ctrl->contig != 1) {
405 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect contig.\n"));
406 return 0;
407 }
408 if (ctrl->minconn != 0 && ctrl->minconn != 1) {
409 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect minconn.\n"));
410 return 0;
411 }
412
413 for (i=0; i<ctrl->ncon; i++) {
414 sum = rsum(ctrl->nparts, ctrl->tpwgts+i, ctrl->ncon);
415 if (sum < 0.99 || sum > 1.01) {
416 IFSET(dbglvl, METIS_DBG_INFO,
417 printf("Input Error: Incorrect sum of %"PRREAL" for tpwgts for constraint %"PRIDX".\n", sum, i));
418 return 0;
419 }
420 }
421 for (i=0; i<ctrl->ncon; i++) {
422 for (j=0; j<ctrl->nparts; j++) {
423 if (ctrl->tpwgts[j*ctrl->ncon+i] <= 0.0) {
424 IFSET(dbglvl, METIS_DBG_INFO,
425 printf("Input Error: Incorrect tpwgts for partition %"PRIDX" and constraint %"PRIDX".\n", j, i));
426 return 0;
427 }
428 }
429 }
430
431 for (i=0; i<ctrl->ncon; i++) {
432 if (ctrl->ubfactors[i] <= 1.0) {
433 IFSET(dbglvl, METIS_DBG_INFO,
434 printf("Input Error: Incorrect ubfactor for constraint %"PRIDX".\n", i));
435 return 0;
436 }
437 }
438
439 break;
440
441
442
443 case METIS_OP_OMETIS:
444 if (ctrl->objtype != METIS_OBJTYPE_NODE) {
445 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect objective type.\n"));
446 return 0;
447 }
448 if (ctrl->ctype != METIS_CTYPE_RM && ctrl->ctype != METIS_CTYPE_SHEM) {
449 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect coarsening scheme.\n"));
450 return 0;
451 }
452 if (ctrl->iptype != METIS_IPTYPE_EDGE && ctrl->iptype != METIS_IPTYPE_NODE) {
453 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect initial partitioning scheme.\n"));
454 return 0;
455 }
456 if (ctrl->rtype != METIS_RTYPE_SEP1SIDED && ctrl->rtype != METIS_RTYPE_SEP2SIDED) {
457 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect refinement scheme.\n"));
458 return 0;
459 }
460 if (ctrl->nseps <= 0) {
461 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect nseps.\n"));
462 return 0;
463 }
464 if (ctrl->niter <= 0) {
465 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect niter.\n"));
466 return 0;
467 }
468 if (ctrl->ufactor <= 0) {
469 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect ufactor.\n"));
470 return 0;
471 }
472 if (ctrl->numflag != 0 && ctrl->numflag != 1) {
473 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect numflag.\n"));
474 return 0;
475 }
476 if (ctrl->nparts != 3) {
477 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect nparts.\n"));
478 return 0;
479 }
480 if (ctrl->ncon != 1) {
481 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect ncon.\n"));
482 return 0;
483 }
484 if (ctrl->compress != 0 && ctrl->compress != 1) {
485 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect compress.\n"));
486 return 0;
487 }
488 if (ctrl->ccorder != 0 && ctrl->ccorder != 1) {
489 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect ccorder.\n"));
490 return 0;
491 }
492 if (ctrl->pfactor < 0.0 ) {
493 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect pfactor.\n"));
494 return 0;
495 }
496
497 for (i=0; i<ctrl->ncon; i++) {
498 if (ctrl->ubfactors[i] <= 1.0) {
499 IFSET(dbglvl, METIS_DBG_INFO,
500 printf("Input Error: Incorrect ubfactor for constraint %"PRIDX".\n", i));
501 return 0;
502 }
503 }
504
505 break;
506
507 default:
508 IFSET(dbglvl, METIS_DBG_INFO, printf("Input Error: Incorrect optype\n"));
509 return 0;
510 }
511
512 return 1;
513 }
514
515
516 /*************************************************************************/
517 /*! This function frees the memory associated with a ctrl_t */
518 /*************************************************************************/
FreeCtrl(ctrl_t ** r_ctrl)519 void FreeCtrl(ctrl_t **r_ctrl)
520 {
521 ctrl_t *ctrl = *r_ctrl;
522
523 FreeWorkSpace(ctrl);
524
525 gk_free((void **)&ctrl->tpwgts, &ctrl->pijbm,
526 &ctrl->ubfactors, &ctrl->maxvwgt, &ctrl, LTERM);
527
528 *r_ctrl = NULL;
529 }
530
531
532