Actual source code: partitioner.c

petsc-master 2020-08-25
Report Typos and Errors
  1:  #include <petsc/private/partitionerimpl.h>

  3: /*@C
  4:   PetscPartitionerSetType - Builds a particular PetscPartitioner

  6:   Collective on PetscPartitioner

  8:   Input Parameters:
  9: + part - The PetscPartitioner object
 10: - name - The kind of partitioner

 12:   Options Database Key:
 13: . -petscpartitioner_type <type> - Sets the PetscPartitioner type; use -help for a list of available types

 15:   Level: intermediate

 17: .seealso: PetscPartitionerGetType(), PetscPartitionerCreate()
 18: @*/
 19: PetscErrorCode PetscPartitionerSetType(PetscPartitioner part, PetscPartitionerType name)
 20: {
 21:   PetscErrorCode (*r)(PetscPartitioner);
 22:   PetscBool      match;

 27:   PetscObjectTypeCompare((PetscObject) part, name, &match);
 28:   if (match) return(0);

 30:   PetscPartitionerRegisterAll();
 31:   PetscFunctionListFind(PetscPartitionerList, name, &r);
 32:   if (!r) SETERRQ1(PetscObjectComm((PetscObject) part), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscPartitioner type: %s", name);

 34:   if (part->ops->destroy) {
 35:     (*part->ops->destroy)(part);
 36:   }
 37:   part->noGraph = PETSC_FALSE;
 38:   PetscMemzero(part->ops, sizeof(*part->ops));
 39:   PetscObjectChangeTypeName((PetscObject) part, name);
 40:   (*r)(part);
 41:   return(0);
 42: }

 44: /*@C
 45:   PetscPartitionerGetType - Gets the PetscPartitioner type name (as a string) from the object.

 47:   Not Collective

 49:   Input Parameter:
 50: . part - The PetscPartitioner

 52:   Output Parameter:
 53: . name - The PetscPartitioner type name

 55:   Level: intermediate

 57: .seealso: PetscPartitionerSetType(), PetscPartitionerCreate()
 58: @*/
 59: PetscErrorCode PetscPartitionerGetType(PetscPartitioner part, PetscPartitionerType *name)
 60: {
 64:   *name = ((PetscObject) part)->type_name;
 65:   return(0);
 66: }

 68: /*@C
 69:    PetscPartitionerViewFromOptions - View from Options

 71:    Collective on PetscPartitioner

 73:    Input Parameters:
 74: +  A - the PetscPartitioner object
 75: .  obj - Optional object
 76: -  name - command line option

 78:    Level: intermediate
 79: .seealso:  PetscPartitionerView(), PetscObjectViewFromOptions()
 80: @*/
 81: PetscErrorCode PetscPartitionerViewFromOptions(PetscPartitioner A,PetscObject obj,const char name[])
 82: {

 87:   PetscObjectViewFromOptions((PetscObject)A,obj,name);
 88:   return(0);
 89: }

 91: /*@
 92:   PetscPartitionerView - Views a PetscPartitioner

 94:   Collective on PetscPartitioner

 96:   Input Parameter:
 97: + part - the PetscPartitioner object to view
 98: - v    - the viewer

100:   Level: developer

102: .seealso: PetscPartitionerDestroy()
103: @*/
104: PetscErrorCode PetscPartitionerView(PetscPartitioner part, PetscViewer v)
105: {
106:   PetscMPIInt    size;
107:   PetscBool      isascii;

112:   if (!v) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) part), &v);}
113:   PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &isascii);
114:   if (isascii) {
115:     MPI_Comm_size(PetscObjectComm((PetscObject) part), &size);
116:     PetscViewerASCIIPrintf(v, "Graph Partitioner: %d MPI Process%s\n", size, size > 1 ? "es" : "");
117:     PetscViewerASCIIPrintf(v, "  type: %s\n", ((PetscObject)part)->type_name);
118:     PetscViewerASCIIPrintf(v, "  edge cut: %D\n", part->edgeCut);
119:     PetscViewerASCIIPrintf(v, "  balance: %.2g\n", part->balance);
120:     PetscViewerASCIIPrintf(v, "  use vertex weights: %d\n", part->usevwgt);
121:   }
122:   if (part->ops->view) {(*part->ops->view)(part, v);}
123:   return(0);
124: }

126: static PetscErrorCode PetscPartitionerGetDefaultType(MPI_Comm comm, const char *currentType, const char **defaultType)
127: {
128:   PetscMPIInt    size;

132:   MPI_Comm_size(comm, &size);
133:   if (!currentType) {
134:     if (size == 1) {
135:       *defaultType = PETSCPARTITIONERSIMPLE;
136:     } else {
137: #if defined(PETSC_HAVE_PARMETIS)
138:       *defaultType = PETSCPARTITIONERPARMETIS;
139: #elif defined(PETSC_HAVE_PTSCOTCH)
140:       *defaultType = PETSCPARTITIONERPTSCOTCH;
141: #elif defined(PETSC_HAVE_CHACO)
142:       *defaultType = PETSCPARTITIONERCHACO;
143: #else
144:       *defaultType = PETSCPARTITIONERSIMPLE;
145: #endif
146:     }
147:   } else {
148:     *defaultType = currentType;
149:   }
150:   return(0);
151: }

153: /*@
154:   PetscPartitionerSetFromOptions - sets parameters in a PetscPartitioner from the options database

156:   Collective on PetscPartitioner

158:   Input Parameter:
159: . part - the PetscPartitioner object to set options for

161:   Options Database Keys:
162: +  -petscpartitioner_type <type> - Sets the PetscPartitioner type; use -help for a list of available types
163: .  -petscpartitioner_use_vertex_weights - Uses weights associated with the graph vertices
164: -  -petscpartitioner_view_graph - View the graph each time PetscPartitionerPartition is called. Viewer can be customized, see PetscOptionsGetViewer()

166:   Level: developer

168: .seealso: PetscPartitionerView(), PetscPartitionerSetType(), PetscPartitionerPartition()
169: @*/
170: PetscErrorCode PetscPartitionerSetFromOptions(PetscPartitioner part)
171: {
172:   const char    *defaultType = NULL;
173:   char           name[256];
174:   PetscBool      flg;

179:   PetscPartitionerGetDefaultType(((PetscObject) part)->comm, ((PetscObject) part)->type_name,&defaultType);
180:   PetscObjectOptionsBegin((PetscObject) part);
181:   PetscOptionsFList("-petscpartitioner_type", "Graph partitioner", "PetscPartitionerSetType", PetscPartitionerList, defaultType, name, sizeof(name), &flg);
182:   if (flg) {
183:     PetscPartitionerSetType(part, name);
184:   } else if (!((PetscObject) part)->type_name) {
185:     PetscPartitionerSetType(part, defaultType);
186:   }
187:   PetscOptionsBool("-petscpartitioner_use_vertex_weights","Use vertex weights","",part->usevwgt,&part->usevwgt,NULL);
188:   if (part->ops->setfromoptions) {
189:     (*part->ops->setfromoptions)(PetscOptionsObject,part);
190:   }
191:   PetscViewerDestroy(&part->viewer);
192:   PetscViewerDestroy(&part->viewerGraph);
193:   PetscOptionsGetViewer(((PetscObject) part)->comm, ((PetscObject) part)->options, ((PetscObject) part)->prefix, "-petscpartitioner_view", &part->viewer, NULL, NULL);
194:   PetscOptionsGetViewer(((PetscObject) part)->comm, ((PetscObject) part)->options, ((PetscObject) part)->prefix, "-petscpartitioner_view_graph", &part->viewerGraph, NULL, &part->viewGraph);
195:   /* process any options handlers added with PetscObjectAddOptionsHandler() */
196:   PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject) part);
197:   PetscOptionsEnd();
198:   return(0);
199: }

201: /*@
202:   PetscPartitionerSetUp - Construct data structures for the PetscPartitioner

204:   Collective on PetscPartitioner

206:   Input Parameter:
207: . part - the PetscPartitioner object to setup

209:   Level: developer

211: .seealso: PetscPartitionerView(), PetscPartitionerDestroy()
212: @*/
213: PetscErrorCode PetscPartitionerSetUp(PetscPartitioner part)
214: {

219:   if (part->ops->setup) {(*part->ops->setup)(part);}
220:   return(0);
221: }

223: /*@
224:   PetscPartitionerReset - Resets data structures for the PetscPartitioner

226:   Collective on PetscPartitioner

228:   Input Parameter:
229: . part - the PetscPartitioner object to reset

231:   Level: developer

233: .seealso: PetscPartitionerSetUp(), PetscPartitionerDestroy()
234: @*/
235: PetscErrorCode PetscPartitionerReset(PetscPartitioner part)
236: {

241:   if (part->ops->reset) {(*part->ops->reset)(part);}
242:   return(0);
243: }

245: /*@
246:   PetscPartitionerDestroy - Destroys a PetscPartitioner object

248:   Collective on PetscPartitioner

250:   Input Parameter:
251: . part - the PetscPartitioner object to destroy

253:   Level: developer

255: .seealso: PetscPartitionerView()
256: @*/
257: PetscErrorCode PetscPartitionerDestroy(PetscPartitioner *part)
258: {

262:   if (!*part) return(0);

265:   if (--((PetscObject)(*part))->refct > 0) {*part = NULL; return(0);}
266:   ((PetscObject) (*part))->refct = 0;

268:   PetscPartitionerReset(*part);

270:   PetscViewerDestroy(&(*part)->viewer);
271:   PetscViewerDestroy(&(*part)->viewerGraph);
272:   if ((*part)->ops->destroy) {(*(*part)->ops->destroy)(*part);}
273:   PetscHeaderDestroy(part);
274:   return(0);
275: }

277: /*@
278:   PetscPartitionerPartition - Partition a graph

280:   Collective on PetscPartitioner

282:   Input Parameters:
283: + part    - The PetscPartitioner
284: . nparts  - Number of partitions
285: . numVertices - Number of vertices in the local part of the graph
286: . start - row pointers for the local part of the graph (CSR style)
287: . adjacency - adjacency list (CSR style)
288: . vertexSection - PetscSection describing the absolute weight of each local vertex (can be NULL)
289: - targetSection - PetscSection describing the absolute weight of each partition (can be NULL)

291:   Output Parameters:
292: + partSection     - The PetscSection giving the division of points by partition
293: - partition       - The list of points by partition

295:   Options Database:
296: + -petscpartitioner_view - View the partitioner information
297: - -petscpartitioner_view_graph - View the graph we are partitioning

299:   Notes:
300:     The chart of the vertexSection (if present) must contain [0,numVertices), with the number of dofs in the section specifying the absolute weight for each vertex.
301:     The chart of the targetSection (if present) must contain [0,nparts), with the number of dofs in the section specifying the absolute weight for each partition. This information must be the same across processes, PETSc does not check it.

303:   Level: developer

305: .seealso PetscPartitionerCreate(), PetscSectionCreate(), PetscSectionSetChart(), PetscSectionSetDof()
306: @*/
307: PetscErrorCode PetscPartitionerPartition(PetscPartitioner part, PetscInt nparts, PetscInt numVertices, PetscInt start[], PetscInt adjacency[], PetscSection vertexSection, PetscSection targetSection, PetscSection partSection, IS *partition)
308: {

314:   if (nparts <= 0) SETERRQ(PetscObjectComm((PetscObject) part), PETSC_ERR_ARG_OUTOFRANGE, "Number of parts must be positive");
315:   if (numVertices < 0) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Number of vertices must be non-negative");
316:   if (numVertices && !part->noGraph) {
320:   }
321:   if (vertexSection) {
322:     PetscInt s,e;

325:     PetscSectionGetChart(vertexSection, &s, &e);
326:     if (s > 0 || e < numVertices) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid vertexSection chart [%D,%D)",s,e);
327:   }
328:   if (targetSection) {
329:     PetscInt s,e;

332:     PetscSectionGetChart(targetSection, &s, &e);
333:     if (s > 0 || e < nparts) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid targetSection chart [%D,%D)",s,e);
334:   }

338:   PetscSectionReset(partSection);
339:   PetscSectionSetChart(partSection, 0, nparts);
340:   if (nparts == 1) { /* quick */
341:     PetscSectionSetDof(partSection, 0, numVertices);
342:     ISCreateStride(PetscObjectComm((PetscObject)part),numVertices,0,1,partition);
343:   } else {
344:     if (!part->ops->partition) SETERRQ1(PetscObjectComm((PetscObject) part), PETSC_ERR_SUP, "PetscPartitioner %s has no partitioning method", ((PetscObject)part)->type_name);
345:     (*part->ops->partition)(part, nparts, numVertices, start, adjacency, vertexSection, targetSection, partSection, partition);
346:   }
347:   PetscSectionSetUp(partSection);
348:   if (part->viewerGraph) {
349:     PetscViewer viewer = part->viewerGraph;
350:     PetscBool   isascii;
351:     PetscInt    v, i;
352:     PetscMPIInt rank;

354:     MPI_Comm_rank(PetscObjectComm((PetscObject) viewer), &rank);
355:     PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &isascii);
356:     if (isascii) {
357:       PetscViewerASCIIPushSynchronized(viewer);
358:       PetscViewerASCIISynchronizedPrintf(viewer, "[%d]Nv: %D\n", rank, numVertices);
359:       for (v = 0; v < numVertices; ++v) {
360:         const PetscInt s = start[v];
361:         const PetscInt e = start[v+1];

363:         PetscViewerASCIISynchronizedPrintf(viewer, "[%d]  ", rank);
364:         for (i = s; i < e; ++i) {PetscViewerASCIISynchronizedPrintf(viewer, "%D ", adjacency[i]);}
365:         PetscViewerASCIISynchronizedPrintf(viewer, "[%D-%D)\n", s, e);
366:       }
367:       PetscViewerFlush(viewer);
368:       PetscViewerASCIIPopSynchronized(viewer);
369:     }
370:   }
371:   if (part->viewer) {
372:     PetscPartitionerView(part,part->viewer);
373:   }
374:   return(0);
375: }

377: /*@
378:   PetscPartitionerCreate - Creates an empty PetscPartitioner object. The type can then be set with PetscPartitionerSetType().

380:   Collective

382:   Input Parameter:
383: . comm - The communicator for the PetscPartitioner object

385:   Output Parameter:
386: . part - The PetscPartitioner object

388:   Level: beginner

390: .seealso: PetscPartitionerSetType(), PETSCPARTITIONERCHACO, PETSCPARTITIONERPARMETIS, PETSCPARTITIONERSHELL, PETSCPARTITIONERSIMPLE, PETSCPARTITIONERGATHER
391: @*/
392: PetscErrorCode PetscPartitionerCreate(MPI_Comm comm, PetscPartitioner *part)
393: {
394:   PetscPartitioner p;
395:   const char       *partitionerType = NULL;
396:   PetscErrorCode   ierr;

400:   *part = NULL;
401:   PetscPartitionerInitializePackage();

403:   PetscHeaderCreate(p, PETSCPARTITIONER_CLASSID, "PetscPartitioner", "Graph Partitioner", "PetscPartitioner", comm, PetscPartitionerDestroy, PetscPartitionerView);
404:   PetscPartitionerGetDefaultType(comm, NULL, &partitionerType);
405:   PetscPartitionerSetType(p,partitionerType);

407:   p->edgeCut = 0;
408:   p->balance = 0.0;
409:   p->usevwgt = PETSC_TRUE;

411:   *part = p;
412:   return(0);
413: }