Actual source code: party.c


  2: #include <../src/mat/impls/adj/mpi/mpiadj.h>

  4: #if defined(PETSC_HAVE_UNISTD_H)
  5:   #include <unistd.h>
  6: #endif

  8: EXTERN_C_BEGIN
  9: #include <party_lib.h>
 10: EXTERN_C_END

 12: typedef struct {
 13:   PetscBool redm;
 14:   PetscBool redo;
 15:   PetscBool recursive;
 16:   PetscBool verbose;
 17:   char      global[15];   /* global method */
 18:   char      local[15];    /* local method */
 19:   PetscInt  nbvtxcoarsed; /* number of vertices for the coarse graph */
 20: } MatPartitioning_Party;

 22: #define SIZE_LOG 10000 /* size of buffer for mesg_log */

 24: static PetscErrorCode MatPartitioningApply_Party(MatPartitioning part, IS *partitioning)
 25: {
 26:   int                    perr;
 27:   PetscInt               i, *parttab, *locals, nb_locals, M, N;
 28:   PetscMPIInt            size, rank;
 29:   Mat                    mat = part->adj, matAdj, matSeq, *A;
 30:   Mat_MPIAdj            *adj;
 31:   MatPartitioning_Party *party = (MatPartitioning_Party *)part->data;
 32:   PetscBool              flg;
 33:   IS                     isrow, iscol;
 34:   int                    n, *edge_p, *edge, *vertex_w, p, *part_party, cutsize, redl, rec;
 35:   const char            *redm, *redo;
 36:   char                  *mesg_log;
 37: #if defined(PETSC_HAVE_UNISTD_H)
 38:   int fd_stdout, fd_pipe[2], count, err;
 39: #endif

 42:   MPI_Comm_size(PetscObjectComm((PetscObject)mat), &size);
 43:   MPI_Comm_rank(PetscObjectComm((PetscObject)mat), &rank);
 44:   PetscObjectTypeCompare((PetscObject)mat, MATMPIADJ, &flg);
 45:   if (size > 1) {
 46:     if (flg) {
 47:       MatMPIAdjToSeq(mat, &matSeq);
 48:     } else {
 49:       PetscInfo(part, "Converting distributed matrix to sequential: this could be a performance loss\n");
 50:       MatGetSize(mat, &M, &N);
 51:       ISCreateStride(PETSC_COMM_SELF, M, 0, 1, &isrow);
 52:       ISCreateStride(PETSC_COMM_SELF, N, 0, 1, &iscol);
 53:       MatCreateSubMatrices(mat, 1, &isrow, &iscol, MAT_INITIAL_MATRIX, &A);
 54:       ISDestroy(&isrow);
 55:       ISDestroy(&iscol);
 56:       matSeq = *A;
 57:       PetscFree(A);
 58:     }
 59:   } else {
 60:     PetscObjectReference((PetscObject)mat);
 61:     matSeq = mat;
 62:   }

 64:   if (!flg) { /* convert regular matrix to MPIADJ */
 65:     MatConvert(matSeq, MATMPIADJ, MAT_INITIAL_MATRIX, &matAdj);
 66:   } else {
 67:     PetscObjectReference((PetscObject)matSeq);
 68:     matAdj = matSeq;
 69:   }

 71:   adj = (Mat_MPIAdj *)matAdj->data; /* finaly adj contains adjacency graph */

 73:   /* arguments for Party library */
 74:   n        = mat->rmap->N;             /* number of vertices in full graph */
 75:   edge_p   = adj->i;                   /* start of edge list for each vertex */
 76:   edge     = adj->j;                   /* edge list data */
 77:   vertex_w = part->vertex_weights;     /* weights for all vertices */
 78:   p        = part->n;                  /* number of parts to create */
 79:   redl     = party->nbvtxcoarsed;      /* how many vertices to coarsen down to? */
 80:   rec      = party->recursive ? 1 : 0; /* recursive bisection */
 81:   redm     = party->redm ? "lam" : ""; /* matching method */
 82:   redo     = party->redo ? "w3" : "";  /* matching optimization method */

 84:   PetscMalloc1(mat->rmap->N, &part_party);

 86:   /* redirect output to buffer */
 87: #if defined(PETSC_HAVE_UNISTD_H)
 88:   fd_stdout = dup(1);
 90:   close(1);
 91:   dup2(fd_pipe[1], 1);
 92:   PetscMalloc1(SIZE_LOG, &mesg_log);
 93: #endif

 95:   /* library call */
 96:   party_lib_times_start();
 97:   perr = party_lib(n, vertex_w, NULL, NULL, NULL, edge_p, edge, NULL, p, part_party, &cutsize, redl, (char *)redm, (char *)redo, party->global, party->local, rec, 1);

 99:   party_lib_times_output(1);
100:   part_info(n, vertex_w, edge_p, edge, NULL, p, part_party, 1);

102: #if defined(PETSC_HAVE_UNISTD_H)
103:   err = fflush(stdout);
105:   count = read(fd_pipe[0], mesg_log, (SIZE_LOG - 1) * sizeof(char));
106:   if (count < 0) count = 0;
107:   mesg_log[count] = 0;
108:   close(1);
109:   dup2(fd_stdout, 1);
110:   close(fd_stdout);
111:   close(fd_pipe[0]);
112:   close(fd_pipe[1]);
113:   if (party->verbose) PetscPrintf(PetscObjectComm((PetscObject)mat), "%s", mesg_log);
114:   PetscFree(mesg_log);
115: #endif

118:   PetscMalloc1(mat->rmap->N, &parttab);
119:   for (i = 0; i < mat->rmap->N; i++) parttab[i] = part_party[i];

121:   /* creation of the index set */
122:   nb_locals = mat->rmap->n;
123:   locals    = parttab + mat->rmap->rstart;

125:   ISCreateGeneral(PetscObjectComm((PetscObject)part), nb_locals, locals, PETSC_COPY_VALUES, partitioning);

127:   /* clean up */
128:   PetscFree(parttab);
129:   PetscFree(part_party);
130:   MatDestroy(&matSeq);
131:   MatDestroy(&matAdj);
132:   return 0;
133: }

135: PetscErrorCode MatPartitioningView_Party(MatPartitioning part, PetscViewer viewer)
136: {
137:   MatPartitioning_Party *party = (MatPartitioning_Party *)part->data;
138:   PetscBool              isascii;

140:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii);
141:   if (isascii) {
142:     PetscViewerASCIIPrintf(viewer, "  Global method: %s\n", party->global);
143:     PetscViewerASCIIPrintf(viewer, "  Local method: %s\n", party->local);
144:     PetscViewerASCIIPrintf(viewer, "  Number of vertices for the coarse graph: %d\n", party->nbvtxcoarsed);
145:     if (party->redm) PetscViewerASCIIPrintf(viewer, "  Using matching method for graph reduction\n");
146:     if (party->redo) PetscViewerASCIIPrintf(viewer, "  Using matching optimization\n");
147:     if (party->recursive) PetscViewerASCIIPrintf(viewer, "  Using recursive bipartitioning\n");
148:   }
149:   return 0;
150: }

152: /*@C
153:    MatPartitioningPartySetGlobal - Set global method for Party partitioner.

155:    Collective

157:    Input Parameters:
158: +  part - the partitioning context
159: -  method - a string representing the method

161:    Options Database Key:
162: .  -mat_partitioning_party_global <method> - the global method

164:    Level: advanced

166:    Note:
167:    The method may be one of `MP_PARTY_OPT`, `MP_PARTY_LIN`, `MP_PARTY_SCA`,
168:    `MP_PARTY_RAN`, `MP_PARTY_GBF`, `MP_PARTY_GCF`, `MP_PARTY_BUB` or `MP_PARTY_DEF`, or
169:    alternatively a string describing the method. Two or more methods can be
170:    combined like "gbf,gcf". Check the Party Library Users Manual for details.

172: .seealso: `MATPARTITIONINGPARTY`, `MatPartitioningPartySetLocal()`
173: @*/
174: PetscErrorCode MatPartitioningPartySetGlobal(MatPartitioning part, const char *global)
175: {
177:   PetscTryMethod(part, "MatPartitioningPartySetGlobal_C", (MatPartitioning, const char *), (part, global));
178:   return 0;
179: }

181: PetscErrorCode MatPartitioningPartySetGlobal_Party(MatPartitioning part, const char *global)
182: {
183:   MatPartitioning_Party *party = (MatPartitioning_Party *)part->data;

185:   PetscStrncpy(party->global, global, 15);
186:   return 0;
187: }

189: /*@C
190:    MatPartitioningPartySetLocal - Set local method used by the Party partitioner.

192:    Collective

194:    Input Parameters:
195: +  part - the partitioning context
196: -  method - a string representing the method

198:    Options Database Key:
199: .  -mat_partitioning_party_local <method> - the local method

201:    Level: advanced

203:    Note:
204:    The method may be one of `MP_PARTY_HELPFUL_SETS`, `MP_PARTY_KERNIGHAN_LIN`, or
205:    `MP_PARTY_NONE`. Check the Party Library Users Manual for details.

207: .seealso: `MATPARTITIONINGPARTY`, `MatPartitioningPartySetGlobal()`
208: @*/
209: PetscErrorCode MatPartitioningPartySetLocal(MatPartitioning part, const char *local)
210: {
212:   PetscTryMethod(part, "MatPartitioningPartySetLocal_C", (MatPartitioning, const char *), (part, local));
213:   return 0;
214: }

216: PetscErrorCode MatPartitioningPartySetLocal_Party(MatPartitioning part, const char *local)

218: {
219:   MatPartitioning_Party *party = (MatPartitioning_Party *)part->data;

221:   PetscStrncpy(party->local, local, 15);
222:   return 0;
223: }

225: /*@
226:    MatPartitioningPartySetCoarseLevel - Set the coarse level parameter for the
227:    Party partitioner.

229:    Collective

231:    Input Parameters:
232: +  part - the partitioning context
233: -  level - the coarse level in range [0.0,1.0]

235:    Options Database Key:
236: .  -mat_partitioning_party_coarse <l> - Coarse level

238:    Level: advanced

240: .seealso: `MATPARTITIONINGPARTY`
241: @*/
242: PetscErrorCode MatPartitioningPartySetCoarseLevel(MatPartitioning part, PetscReal level)
243: {
246:   PetscTryMethod(part, "MatPartitioningPartySetCoarseLevel_C", (MatPartitioning, PetscReal), (part, level));
247:   return 0;
248: }

250: PetscErrorCode MatPartitioningPartySetCoarseLevel_Party(MatPartitioning part, PetscReal level)
251: {
252:   MatPartitioning_Party *party = (MatPartitioning_Party *)part->data;

255:   party->nbvtxcoarsed = (PetscInt)(part->adj->cmap->N * level);
256:   if (party->nbvtxcoarsed < 20) party->nbvtxcoarsed = 20;
257:   return 0;
258: }

260: /*@
261:    MatPartitioningPartySetMatchOptimization - Activate matching optimization for
262:    graph reduction.

264:    Collective

266:    Input Parameters:
267: +  part - the partitioning context
268: -  opt - boolean flag

270:    Options Database Key:
271: .  -mat_partitioning_party_match_optimization - Matching optimization on/off

273:    Level: advanced

275: .seealso:  `MATPARTITIONINGPARTY`
276: @*/
277: PetscErrorCode MatPartitioningPartySetMatchOptimization(MatPartitioning part, PetscBool opt)
278: {
281:   PetscTryMethod(part, "MatPartitioningPartySetMatchOptimization_C", (MatPartitioning, PetscBool), (part, opt));
282:   return 0;
283: }

285: PetscErrorCode MatPartitioningPartySetMatchOptimization_Party(MatPartitioning part, PetscBool opt)
286: {
287:   MatPartitioning_Party *party = (MatPartitioning_Party *)part->data;

289:   party->redo = opt;
290:   return 0;
291: }

293: /*@
294:    MatPartitioningPartySetBipart - Activate or deactivate recursive bisection in the Party partitioner

296:    Collective

298:    Input Parameters:
299: +  part - the partitioning context
300: -  bp - boolean flag

302:    Options Database Key:
303: -  -mat_partitioning_party_bipart - Bipartitioning option on/off

305:    Level: advanced

307: .seealso:  `MATPARTITIONINGPARTY`
308: @*/
309: PetscErrorCode MatPartitioningPartySetBipart(MatPartitioning part, PetscBool bp)
310: {
313:   PetscTryMethod(part, "MatPartitioningPartySetBipart_C", (MatPartitioning, PetscBool), (part, bp));
314:   return 0;
315: }

317: PetscErrorCode MatPartitioningPartySetBipart_Party(MatPartitioning part, PetscBool bp)
318: {
319:   MatPartitioning_Party *party = (MatPartitioning_Party *)part->data;

321:   party->recursive = bp;
322:   return 0;
323: }

325: PetscErrorCode MatPartitioningSetFromOptions_Party(MatPartitioning part, PetscOptionItems *PetscOptionsObject)
326: {
327:   PetscBool              flag;
328:   char                   value[256];
329:   PetscReal              r;
330:   MatPartitioning_Party *party = (MatPartitioning_Party *)part->data;

332:   PetscOptionsHeadBegin(PetscOptionsObject, "Set Party partitioning options");
333:   PetscOptionsString("-mat_partitioning_party_global", "Global method", "MatPartitioningPartySetGlobal", party->global, value, sizeof(value), &flag);
334:   if (flag) MatPartitioningPartySetGlobal(part, value);
335:   PetscOptionsString("-mat_partitioning_party_local", "Local method", "MatPartitioningPartySetLocal", party->local, value, sizeof(value), &flag);
336:   if (flag) MatPartitioningPartySetLocal(part, value);
337:   PetscOptionsReal("-mat_partitioning_party_coarse", "Coarse level", "MatPartitioningPartySetCoarseLevel", 0.0, &r, &flag);
338:   if (flag) MatPartitioningPartySetCoarseLevel(part, r);
339:   PetscOptionsBool("-mat_partitioning_party_match_optimization", "Matching optimization on/off", "MatPartitioningPartySetMatchOptimization", party->redo, &party->redo, NULL);
340:   PetscOptionsBool("-mat_partitioning_party_bipart", "Bipartitioning on/off", "MatPartitioningPartySetBipart", party->recursive, &party->recursive, NULL);
341:   PetscOptionsBool("-mat_partitioning_party_verbose", "Show library output", "", party->verbose, &party->verbose, NULL);
342:   PetscOptionsHeadEnd();
343:   return 0;
344: }

346: PetscErrorCode MatPartitioningDestroy_Party(MatPartitioning part)
347: {
348:   MatPartitioning_Party *party = (MatPartitioning_Party *)part->data;

350:   PetscFree(party);
351:   /* clear composed functions */
352:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetGlobal_C", NULL);
353:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetLocal_C", NULL);
354:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetCoarseLevel_C", NULL);
355:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetMatchOptimization_C", NULL);
356:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetBipart_C", NULL);
357:   return 0;
358: }

360: /*MC
361:    MATPARTITIONINGPARTY - Creates a partitioning context via the external package Party.

363:    Level: beginner

365:    Notes:
366:     See http://wwwcs.upb.de/fachbereich/AG/monien/RESEARCH/PART/party.html

368:     Does not support the `MatPartitioningSetUseEdgeWeights()` option

370: .seealso: `MatPartitioningSetType()`, `MatPartitioningType`, `MatPartitioningPartySetGlobal()`, `MatPartitioningPartySetLocal()`,
371:           `MatPartitioningPartySetCoarseLevel()`, `MatPartitioningPartySetMatchOptimization()`, `MatPartitioningPartySetBipart()`
372: M*/

374: PETSC_EXTERN PetscErrorCode MatPartitioningCreate_Party(MatPartitioning part)
375: {
376:   MatPartitioning_Party *party;

378:   PetscNew(&party);
379:   part->data = (void *)party;

381:   PetscStrcpy(party->global, "gcf,gbf");
382:   PetscStrcpy(party->local, "kl");

384:   party->redm         = PETSC_TRUE;
385:   party->redo         = PETSC_TRUE;
386:   party->recursive    = PETSC_TRUE;
387:   party->verbose      = PETSC_FALSE;
388:   party->nbvtxcoarsed = 200;

390:   part->ops->apply          = MatPartitioningApply_Party;
391:   part->ops->view           = MatPartitioningView_Party;
392:   part->ops->destroy        = MatPartitioningDestroy_Party;
393:   part->ops->setfromoptions = MatPartitioningSetFromOptions_Party;

395:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetGlobal_C", MatPartitioningPartySetGlobal_Party);
396:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetLocal_C", MatPartitioningPartySetLocal_Party);
397:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetCoarseLevel_C", MatPartitioningPartySetCoarseLevel_Party);
398:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetMatchOptimization_C", MatPartitioningPartySetMatchOptimization_Party);
399:   PetscObjectComposeFunction((PetscObject)part, "MatPartitioningPartySetBipart_C", MatPartitioningPartySetBipart_Party);
400:   return 0;
401: }