1 /*
2 * Copyright 1997, Regents of the University of Minnesota
3 *
4 * main.c
5 *
6 * This file contains code for testing teh adaptive partitioning routines
7 *
8 * Started 5/19/97
9 * George
10 *
11 * $Id: ptest.c,v 1.3 2003/07/22 21:47:20 karypis Exp $
12 *
13 */
14
15 #include <parmetisbin.h>
16
17
18 /*************************************************************************/
19 /*! Entry point of the testing routine */
20 /*************************************************************************/
main(idx_t argc,char * argv[])21 idx_t main(idx_t argc, char *argv[])
22 {
23 idx_t mype, npes;
24 MPI_Comm comm;
25
26 MPI_Init(&argc, &argv);
27 MPI_Comm_dup(MPI_COMM_WORLD, &comm);
28 MPI_Comm_size(comm, &npes);
29 MPI_Comm_rank(comm, &mype);
30
31 if (argc != 2) {
32 if (mype == 0)
33 printf("Usage: %s <graph-file>\n", argv[0]);
34
35 MPI_Finalize();
36 exit(0);
37 }
38
39 TestParMetis_GPart(argv[1], comm);
40
41 MPI_Comm_free(&comm);
42
43 MPI_Finalize();
44
45 return 0;
46 }
47
48
49
50 /***********************************************************************************/
51 /*! This function tests the various graph partitioning and ordering routines */
52 /***********************************************************************************/
TestParMetis(char * filename,MPI_Comm comm)53 void TestParMetis(char *filename, MPI_Comm comm)
54 {
55 idx_t ncon, nparts, npes, mype, opt2, realcut;
56 graph_t graph, mgraph;
57 idx_t *part, *mpart, *savepart, *order, *sizes;
58 idx_t numflag=0, wgtflag=0, options[10], edgecut, ndims;
59 real_t ipc2redist, *xyz, *tpwgts = NULL, ubvec[MAXNCON];
60
61 MPI_Comm_size(comm, &npes);
62 MPI_Comm_rank(comm, &mype);
63
64 ndims = 2;
65
66 ParallelReadGraph(&graph, filename, comm);
67 xyz = ReadTestCoordinates(&graph, filename, 2, comm);
68 MPI_Barrier(comm);
69
70 part = imalloc(graph.nvtxs, "TestParMetis_V3: part");
71 tpwgts = rmalloc(MAXNCON*npes*2, "TestParMetis_V3: tpwgts");
72 rset(MAXNCON, 1.05, ubvec);
73 graph.vwgt = ismalloc(graph.nvtxs*5, 1, "TestParMetis_V3: vwgt");
74
75
76 /*======================================================================
77 / ParMETIS_V3_PartKway
78 /=======================================================================*/
79 options[0] = 1;
80 options[1] = 3;
81 options[2] = 1;
82 wgtflag = 2;
83 numflag = 0;
84 edgecut = 0;
85 for (nparts=2*npes; nparts>=npes/2 && nparts > 0; nparts = nparts/2) {
86 for (ncon=1; ncon<=5; ncon+=2) {
87
88 if (ncon > 1 && nparts > 1)
89 Mc_AdaptGraph(&graph, part, ncon, nparts, comm);
90 else
91 iset(graph.nvtxs, 1, graph.vwgt);
92
93 for (opt2=1; opt2<=2; opt2++) {
94 options[2] = opt2;
95
96 rset(nparts*ncon, 1.0/(real_t)nparts, tpwgts);
97 if (mype == 0)
98 printf("\nTesting ParMETIS_V3_PartKway with options[1-2] = {%"PRIDX" %"PRIDX"}, Ncon: %"PRIDX", Nparts: %"PRIDX"\n", options[1], options[2], ncon, nparts);
99
100 ParMETIS_V3_PartKway(graph.vtxdist, graph.xadj, graph.adjncy, graph.vwgt, NULL, &wgtflag,
101 &numflag, &ncon, &nparts, tpwgts, ubvec, options, &edgecut, part, &comm);
102
103 if (mype == 0) {
104 printf("ParMETIS_V3_PartKway reported a cut of %"PRIDX"\n", edgecut);
105 }
106 }
107 }
108 }
109
110
111 /*======================================================================
112 / ParMETIS_V3_PartGeomKway
113 /=======================================================================*/
114 options[0] = 1;
115 options[1] = 3;
116 wgtflag = 2;
117 numflag = 0;
118 for (nparts=2*npes; nparts>=npes/2 && nparts > 0; nparts = nparts/2) {
119 for (ncon=1; ncon<=5; ncon+=2) {
120
121 if (ncon > 1)
122 Mc_AdaptGraph(&graph, part, ncon, nparts, comm);
123 else
124 iset(graph.nvtxs, 1, graph.vwgt);
125
126 for (opt2=1; opt2<=2; opt2++) {
127 options[2] = opt2;
128
129 rset(nparts*ncon, 1.0/(real_t)nparts, tpwgts);
130 if (mype == 0)
131 printf("\nTesting ParMETIS_V3_PartGeomKway with options[1-2] = {%"PRIDX" %"PRIDX"}, Ncon: %"PRIDX", Nparts: %"PRIDX"\n", options[1], options[2], ncon, nparts);
132
133 ParMETIS_V3_PartGeomKway(graph.vtxdist, graph.xadj, graph.adjncy, graph.vwgt, NULL, &wgtflag,
134 &numflag, &ndims, xyz, &ncon, &nparts, tpwgts, ubvec, options, &edgecut, part, &comm);
135
136 if (mype == 0) {
137 printf("ParMETIS_V3_PartGeomKway reported a cut of %"PRIDX"\n", edgecut);
138 }
139 }
140 }
141 }
142
143
144
145 /*======================================================================
146 / ParMETIS_V3_PartGeom
147 /=======================================================================*/
148 wgtflag = 0;
149 numflag = 0;
150 if (mype == 0)
151 printf("\nTesting ParMETIS_V3_PartGeom\n");
152
153 /* ParMETIS_V3_PartGeom(graph.vtxdist, &ndims, xyz, part, &comm); */
154
155 if (mype == 0)
156 printf("ParMETIS_V3_PartGeom partition complete\n");
157 /*
158 realcut = ComputeRealCut(graph.vtxdist, part, filename, comm);
159 if (mype == 0)
160 printf("ParMETIS_V3_PartGeom reported a cut of %"PRIDX"\n", realcut);
161 */
162
163 /*======================================================================
164 / ParMETIS_V3_RefineKway
165 /=======================================================================*/
166 options[0] = 1;
167 options[1] = 3;
168 options[2] = 1;
169 options[3] = COUPLED;
170 nparts = npes;
171 wgtflag = 0;
172 numflag = 0;
173 ncon = 1;
174 rset(nparts*ncon, 1.0/(real_t)nparts, tpwgts);
175
176 if (mype == 0)
177 printf("\nTesting ParMETIS_V3_RefineKway with default options (before move)\n");
178
179 ParMETIS_V3_RefineKway(graph.vtxdist, graph.xadj, graph.adjncy, NULL, NULL, &wgtflag,
180 &numflag, &ncon, &nparts, tpwgts, ubvec, options, &edgecut, part, &comm);
181
182 MALLOC_CHECK(NULL);
183
184 if (mype == 0) {
185 printf("ParMETIS_V3_RefineKway reported a cut of %"PRIDX"\n", edgecut);
186 }
187
188
189 MALLOC_CHECK(NULL);
190
191 /* Compute a good partition and move the graph. Do so quietly! */
192 options[0] = 0;
193 nparts = npes;
194 wgtflag = 0;
195 numflag = 0;
196 ncon = 1;
197 rset(nparts*ncon, 1.0/(real_t)nparts, tpwgts);
198 ParMETIS_V3_PartKway(graph.vtxdist, graph.xadj, graph.adjncy, NULL, NULL, &wgtflag,
199 &numflag, &ncon, &npes, tpwgts, ubvec, options, &edgecut, part, &comm);
200 TestMoveGraph(&graph, &mgraph, part, comm);
201 gk_free((void **)&(graph.vwgt), LTERM);
202 mpart = ismalloc(mgraph.nvtxs, mype, "TestParMetis_V3: mpart");
203 savepart = imalloc(mgraph.nvtxs, "TestParMetis_V3: savepart");
204
205 MALLOC_CHECK(NULL);
206
207 /*======================================================================
208 / ParMETIS_V3_RefineKway
209 /=======================================================================*/
210 options[0] = 1;
211 options[1] = 3;
212 options[3] = COUPLED;
213 nparts = npes;
214 wgtflag = 0;
215 numflag = 0;
216
217 for (ncon=1; ncon<=5; ncon+=2) {
218 for (opt2=1; opt2<=2; opt2++) {
219 options[2] = opt2;
220
221 rset(nparts*ncon, 1.0/(real_t)nparts, tpwgts);
222 if (mype == 0)
223 printf("\nTesting ParMETIS_V3_RefineKway with options[1-3] = {%"PRIDX" %"PRIDX" %"PRIDX"}, Ncon: %"PRIDX", Nparts: %"PRIDX"\n", options[1], options[2], options[3], ncon, nparts);
224 ParMETIS_V3_RefineKway(mgraph.vtxdist, mgraph.xadj, mgraph.adjncy, NULL, NULL, &wgtflag,
225 &numflag, &ncon, &nparts, tpwgts, ubvec, options, &edgecut, mpart, &comm);
226
227 if (mype == 0) {
228 printf("ParMETIS_V3_RefineKway reported a cut of %"PRIDX"\n", edgecut);
229 }
230 }
231 }
232
233
234 /*======================================================================
235 / ParMETIS_V3_AdaptiveRepart
236 /=======================================================================*/
237 mgraph.vwgt = ismalloc(mgraph.nvtxs*5, 1, "TestParMetis_V3: mgraph.vwgt");
238 mgraph.vsize = ismalloc(mgraph.nvtxs, 1, "TestParMetis_V3: mgraph.vsize");
239 AdaptGraph(&mgraph, 4, comm);
240 options[0] = 1;
241 options[1] = 7;
242 options[3] = COUPLED;
243 wgtflag = 2;
244 numflag = 0;
245
246 for (nparts=2*npes; nparts>=npes/2; nparts = nparts/2) {
247
248 ncon = 1;
249 wgtflag = 0;
250 options[0] = 0;
251 rset(nparts*ncon, 1.0/(real_t)nparts, tpwgts);
252 ParMETIS_V3_PartKway(mgraph.vtxdist, mgraph.xadj, mgraph.adjncy, NULL, NULL,
253 &wgtflag, &numflag, &ncon, &nparts, tpwgts, ubvec, options, &edgecut, savepart, &comm);
254 options[0] = 1;
255 wgtflag = 2;
256
257 for (ncon=1; ncon<=3; ncon+=2) {
258 rset(nparts*ncon, 1.0/(real_t)nparts, tpwgts);
259
260 if (ncon > 1)
261 Mc_AdaptGraph(&mgraph, savepart, ncon, nparts, comm);
262 else
263 AdaptGraph(&mgraph, 4, comm);
264 /* iset(mgraph.nvtxs, 1, mgraph.vwgt); */
265
266 for (ipc2redist=1000.0; ipc2redist>=0.001; ipc2redist/=1000.0) {
267 for (opt2=1; opt2<=2; opt2++) {
268 icopy(mgraph.nvtxs, savepart, mpart);
269 options[2] = opt2;
270
271 if (mype == 0)
272 printf("\nTesting ParMETIS_V3_AdaptiveRepart with options[1-3] = {%"PRIDX" %"PRIDX" %"PRIDX"}, ipc2redist: %.3"PRREAL", Ncon: %"PRIDX", Nparts: %"PRIDX"\n", options[1], options[2], options[3], ipc2redist, ncon, nparts);
273
274 ParMETIS_V3_AdaptiveRepart(mgraph.vtxdist, mgraph.xadj, mgraph.adjncy, mgraph.vwgt,
275 mgraph.vsize, NULL, &wgtflag, &numflag, &ncon, &nparts, tpwgts, ubvec, &ipc2redist,
276 options, &edgecut, mpart, &comm);
277
278 if (mype == 0) {
279 printf("ParMETIS_V3_AdaptiveRepart reported a cut of %"PRIDX"\n", edgecut);
280 }
281 }
282 }
283 }
284 }
285
286 free(mgraph.vwgt);
287 free(mgraph.vsize);
288
289
290
291 /*======================================================================
292 / ParMETIS_V3_NodeND
293 /=======================================================================*/
294 sizes = imalloc(2*npes, "TestParMetis_V3: sizes");
295 order = imalloc(graph.nvtxs, "TestParMetis_V3: sizes");
296
297 options[0] = 1;
298 options[PMV3_OPTION_DBGLVL] = 3;
299 options[PMV3_OPTION_SEED] = 1;
300 numflag = 0;
301
302 for (opt2=1; opt2<=2; opt2++) {
303 options[PMV3_OPTION_IPART] = opt2;
304
305 if (mype == 0)
306 printf("\nTesting ParMETIS_V3_NodeND with options[1-3] = {%"PRIDX" %"PRIDX" %"PRIDX"}\n", options[1], options[2], options[3]);
307
308 ParMETIS_V3_NodeND(graph.vtxdist, graph.xadj, graph.adjncy, &numflag, options,
309 order, sizes, &comm);
310 }
311
312
313 gk_free((void **)&tpwgts, &part, &mpart, &savepart, &order, &sizes, LTERM);
314
315 }
316
317
318
319 /******************************************************************************
320 * This function takes a partition vector that is distributed and reads in
321 * the original graph and computes the edgecut
322 *******************************************************************************/
ComputeRealCut(idx_t * vtxdist,idx_t * part,char * filename,MPI_Comm comm)323 idx_t ComputeRealCut(idx_t *vtxdist, idx_t *part, char *filename, MPI_Comm comm)
324 {
325 idx_t i, j, nvtxs, mype, npes, cut;
326 idx_t *xadj, *adjncy, *gpart;
327 MPI_Status status;
328
329 MPI_Comm_size(comm, &npes);
330 MPI_Comm_rank(comm, &mype);
331
332 if (mype != 0) {
333 MPI_Send((void *)part, vtxdist[mype+1]-vtxdist[mype], IDX_T, 0, 1, comm);
334 }
335 else { /* Processor 0 does all the rest */
336 gpart = imalloc(vtxdist[npes], "ComputeRealCut: gpart");
337 icopy(vtxdist[1], part, gpart);
338
339 for (i=1; i<npes; i++)
340 MPI_Recv((void *)(gpart+vtxdist[i]), vtxdist[i+1]-vtxdist[i], IDX_T, i, 1, comm, &status);
341
342 ReadMetisGraph(filename, &nvtxs, &xadj, &adjncy);
343
344 /* OK, now compute the cut */
345 for (cut=0, i=0; i<nvtxs; i++) {
346 for (j=xadj[i]; j<xadj[i+1]; j++) {
347 if (gpart[i] != gpart[adjncy[j]])
348 cut++;
349 }
350 }
351 cut = cut/2;
352
353 gk_free((void **)&gpart, &xadj, &adjncy, LTERM);
354
355 return cut;
356 }
357 return 0;
358 }
359
360
361 /******************************************************************************
362 * This function takes a partition vector that is distributed and reads in
363 * the original graph and computes the edgecut
364 *******************************************************************************/
ComputeRealCut2(idx_t * vtxdist,idx_t * mvtxdist,idx_t * part,idx_t * mpart,char * filename,MPI_Comm comm)365 idx_t ComputeRealCut2(idx_t *vtxdist, idx_t *mvtxdist, idx_t *part, idx_t *mpart, char *filename, MPI_Comm comm)
366 {
367 idx_t i, j, nvtxs, mype, npes, cut;
368 idx_t *xadj, *adjncy, *gpart, *gmpart, *perm, *sizes;
369 MPI_Status status;
370
371
372 MPI_Comm_size(comm, &npes);
373 MPI_Comm_rank(comm, &mype);
374
375 if (mype != 0) {
376 MPI_Send((void *)part, vtxdist[mype+1]-vtxdist[mype], IDX_T, 0, 1, comm);
377 MPI_Send((void *)mpart, mvtxdist[mype+1]-mvtxdist[mype], IDX_T, 0, 1, comm);
378 }
379 else { /* Processor 0 does all the rest */
380 gpart = imalloc(vtxdist[npes], "ComputeRealCut: gpart");
381 icopy(vtxdist[1], part, gpart);
382 gmpart = imalloc(mvtxdist[npes], "ComputeRealCut: gmpart");
383 icopy(mvtxdist[1], mpart, gmpart);
384
385 for (i=1; i<npes; i++) {
386 MPI_Recv((void *)(gpart+vtxdist[i]), vtxdist[i+1]-vtxdist[i], IDX_T, i, 1, comm, &status);
387 MPI_Recv((void *)(gmpart+mvtxdist[i]), mvtxdist[i+1]-mvtxdist[i], IDX_T, i, 1, comm, &status);
388 }
389
390 /* OK, now go and reconstruct the permutation to go from the graph to mgraph */
391 perm = imalloc(vtxdist[npes], "ComputeRealCut: perm");
392 sizes = ismalloc(npes+1, 0, "ComputeRealCut: sizes");
393
394 for (i=0; i<vtxdist[npes]; i++)
395 sizes[gpart[i]]++;
396 MAKECSR(i, npes, sizes);
397 for (i=0; i<vtxdist[npes]; i++)
398 perm[i] = sizes[gpart[i]]++;
399
400 /* Ok, now read the graph from the file */
401 ReadMetisGraph(filename, &nvtxs, &xadj, &adjncy);
402
403 /* OK, now compute the cut */
404 for (cut=0, i=0; i<nvtxs; i++) {
405 for (j=xadj[i]; j<xadj[i+1]; j++) {
406 if (gmpart[perm[i]] != gmpart[perm[adjncy[j]]])
407 cut++;
408 }
409 }
410 cut = cut/2;
411
412 gk_free((void **)&gpart, &gmpart, &perm, &sizes, &xadj, &adjncy, LTERM);
413
414 return cut;
415 }
416
417 return 0;
418 }
419
420
421
422 /******************************************************************************
423 * This function takes a graph and its partition vector and creates a new
424 * graph corresponding to the one after the movement
425 *******************************************************************************/
TestMoveGraph(graph_t * ograph,graph_t * omgraph,idx_t * part,MPI_Comm comm)426 void TestMoveGraph(graph_t *ograph, graph_t *omgraph, idx_t *part, MPI_Comm comm)
427 {
428 idx_t npes, mype;
429 ctrl_t ctrl;
430 wspace_t wspace;
431 graph_t *graph, *mgraph;
432 idx_t options[5] = {0, 0, 1, 0, 0};
433
434 MPI_Comm_size(comm, &npes);
435 MPI_Comm_rank(comm, &mype);
436
437 SetUpCtrl(&ctrl, npes, 0, comm);
438 ctrl.CoarsenTo = 1; /* Needed by SetUpGraph, otherwise we can FP errors */
439 graph = SetUpGraph(&ctrl, ograph->vtxdist, ograph->xadj, NULL, ograph->adjncy, NULL, 0);
440 AllocateWSpace(&ctrl, graph, &wspace);
441
442 SetUp(&ctrl, graph, &wspace);
443 graph->where = part;
444 graph->ncon = 1;
445 mgraph = MoveGraph(&ctrl, graph, &wspace);
446
447 omgraph->gnvtxs = mgraph->gnvtxs;
448 omgraph->nvtxs = mgraph->nvtxs;
449 omgraph->nedges = mgraph->nedges;
450 omgraph->vtxdist = mgraph->vtxdist;
451 omgraph->xadj = mgraph->xadj;
452 omgraph->adjncy = mgraph->adjncy;
453 mgraph->vtxdist = NULL;
454 mgraph->xadj = NULL;
455 mgraph->adjncy = NULL;
456 FreeGraph(mgraph);
457
458 graph->where = NULL;
459 FreeInitialGraphAndRemap(graph, 0, 1);
460 FreeWSpace(&wspace);
461 }
462
463
464 /*****************************************************************************
465 * This function sets up a graph data structure for partitioning
466 *****************************************************************************/
SetUpGraph(ctrl_t * ctrl,idx_t * vtxdist,idx_t * xadj,idx_t * vwgt,idx_t * adjncy,idx_t * adjwgt,idx_t wgtflag)467 graph_t *SetUpGraph(ctrl_t *ctrl, idx_t *vtxdist, idx_t *xadj,
468 idx_t *vwgt, idx_t *adjncy, idx_t *adjwgt, idx_t wgtflag)
469 {
470 idx_t mywgtflag;
471
472 mywgtflag = wgtflag;
473 return Mc_SetUpGraph(ctrl, 1, vtxdist, xadj, vwgt, adjncy, adjwgt, &mywgtflag);
474 }
475
476
477