Actual source code: partitioner.c
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: Note:
16: $ PETSCPARTITIONERCHACO - The Chaco partitioner (--download-chaco)
17: $ PETSCPARTITIONERPARMETIS - The ParMetis partitioner (--download-parmetis)
18: $ PETSCPARTITIONERSHELL - A shell partitioner implemented by the user
19: $ PETSCPARTITIONERSIMPLE - A simple partitioner that divides cells into equal, contiguous chunks
20: $ PETSCPARTITIONERGATHER - Gathers all cells onto process 0
22: Level: intermediate
24: .seealso: `PetscPartitionerGetType()`, `PetscPartitionerCreate()`
25: @*/
26: PetscErrorCode PetscPartitionerSetType(PetscPartitioner part, PetscPartitionerType name)
27: {
28: PetscErrorCode (*r)(PetscPartitioner);
29: PetscBool match;
31: PetscFunctionBegin;
33: PetscCall(PetscObjectTypeCompare((PetscObject)part, name, &match));
34: if (match) PetscFunctionReturn(PETSC_SUCCESS);
36: PetscCall(PetscPartitionerRegisterAll());
37: PetscCall(PetscFunctionListFind(PetscPartitionerList, name, &r));
38: PetscCheck(r, PetscObjectComm((PetscObject)part), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscPartitioner type: %s", name);
40: PetscTryTypeMethod(part, destroy);
41: part->noGraph = PETSC_FALSE;
42: PetscCall(PetscMemzero(part->ops, sizeof(*part->ops)));
43: PetscCall(PetscObjectChangeTypeName((PetscObject)part, name));
44: PetscCall((*r)(part));
45: PetscFunctionReturn(PETSC_SUCCESS);
46: }
48: /*@C
49: PetscPartitionerGetType - Gets the PetscPartitioner type name (as a string) from the object.
51: Not Collective
53: Input Parameter:
54: . part - The PetscPartitioner
56: Output Parameter:
57: . name - The PetscPartitioner type name
59: Level: intermediate
61: .seealso: `PetscPartitionerSetType()`, `PetscPartitionerCreate()`
62: @*/
63: PetscErrorCode PetscPartitionerGetType(PetscPartitioner part, PetscPartitionerType *name)
64: {
65: PetscFunctionBegin;
68: *name = ((PetscObject)part)->type_name;
69: PetscFunctionReturn(PETSC_SUCCESS);
70: }
72: /*@C
73: PetscPartitionerViewFromOptions - View from Options
75: Collective on PetscPartitioner
77: Input Parameters:
78: + A - the PetscPartitioner object
79: . obj - Optional object
80: - name - command line option
82: Level: intermediate
83: .seealso: `PetscPartitionerView()`, `PetscObjectViewFromOptions()`
84: @*/
85: PetscErrorCode PetscPartitionerViewFromOptions(PetscPartitioner A, PetscObject obj, const char name[])
86: {
87: PetscFunctionBegin;
89: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
90: PetscFunctionReturn(PETSC_SUCCESS);
91: }
93: /*@
94: PetscPartitionerView - Views a PetscPartitioner
96: Collective on PetscPartitioner
98: Input Parameters:
99: + part - the PetscPartitioner object to view
100: - v - the viewer
102: Level: developer
104: .seealso: `PetscPartitionerDestroy()`
105: @*/
106: PetscErrorCode PetscPartitionerView(PetscPartitioner part, PetscViewer v)
107: {
108: PetscMPIInt size;
109: PetscBool isascii;
111: PetscFunctionBegin;
113: if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)part), &v));
114: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &isascii));
115: if (isascii) {
116: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)part), &size));
117: PetscCall(PetscViewerASCIIPrintf(v, "Graph Partitioner: %d MPI Process%s\n", size, size > 1 ? "es" : ""));
118: PetscCall(PetscViewerASCIIPrintf(v, " type: %s\n", ((PetscObject)part)->type_name));
119: PetscCall(PetscViewerASCIIPrintf(v, " edge cut: %" PetscInt_FMT "\n", part->edgeCut));
120: PetscCall(PetscViewerASCIIPrintf(v, " balance: %.2g\n", (double)part->balance));
121: PetscCall(PetscViewerASCIIPrintf(v, " use vertex weights: %d\n", part->usevwgt));
122: }
123: PetscTryTypeMethod(part, view, v);
124: PetscFunctionReturn(PETSC_SUCCESS);
125: }
127: static PetscErrorCode PetscPartitionerGetDefaultType(MPI_Comm comm, const char **defaultType)
128: {
129: PetscMPIInt size;
131: PetscFunctionBegin;
132: PetscCallMPI(MPI_Comm_size(comm, &size));
133: if (size == 1) {
134: *defaultType = PETSCPARTITIONERSIMPLE;
135: } else {
136: #if defined(PETSC_HAVE_PARMETIS)
137: *defaultType = PETSCPARTITIONERPARMETIS;
138: #elif defined(PETSC_HAVE_PTSCOTCH)
139: *defaultType = PETSCPARTITIONERPTSCOTCH;
140: #elif defined(PETSC_HAVE_CHACO)
141: *defaultType = PETSCPARTITIONERCHACO;
142: #else
143: *defaultType = PETSCPARTITIONERSIMPLE;
144: #endif
145: }
146: PetscFunctionReturn(PETSC_SUCCESS);
147: }
149: /*@
150: PetscPartitionerSetFromOptions - sets parameters in a PetscPartitioner from the options database
152: Collective on PetscPartitioner
154: Input Parameter:
155: . part - the PetscPartitioner object to set options for
157: Options Database Keys:
158: + -petscpartitioner_type <type> - Sets the PetscPartitioner type; use -help for a list of available types
159: . -petscpartitioner_use_vertex_weights - Uses weights associated with the graph vertices
160: - -petscpartitioner_view_graph - View the graph each time PetscPartitionerPartition is called. Viewer can be customized, see PetscOptionsGetViewer()
162: Level: developer
164: .seealso: `PetscPartitionerView()`, `PetscPartitionerSetType()`, `PetscPartitionerPartition()`
165: @*/
166: PetscErrorCode PetscPartitionerSetFromOptions(PetscPartitioner part)
167: {
168: const char *currentType = NULL;
169: char name[256];
170: PetscBool flg;
172: PetscFunctionBegin;
174: PetscObjectOptionsBegin((PetscObject)part);
175: PetscCall(PetscPartitionerGetType(part, ¤tType));
176: PetscCall(PetscOptionsFList("-petscpartitioner_type", "Graph partitioner", "PetscPartitionerSetType", PetscPartitionerList, currentType, name, sizeof(name), &flg));
177: if (flg) PetscCall(PetscPartitionerSetType(part, name));
178: PetscCall(PetscOptionsBool("-petscpartitioner_use_vertex_weights", "Use vertex weights", "", part->usevwgt, &part->usevwgt, NULL));
179: PetscTryTypeMethod(part, setfromoptions, PetscOptionsObject);
180: PetscCall(PetscViewerDestroy(&part->viewer));
181: PetscCall(PetscViewerDestroy(&part->viewerGraph));
182: PetscCall(PetscOptionsGetViewer(((PetscObject)part)->comm, ((PetscObject)part)->options, ((PetscObject)part)->prefix, "-petscpartitioner_view", &part->viewer, NULL, NULL));
183: PetscCall(PetscOptionsGetViewer(((PetscObject)part)->comm, ((PetscObject)part)->options, ((PetscObject)part)->prefix, "-petscpartitioner_view_graph", &part->viewerGraph, NULL, &part->viewGraph));
184: /* process any options handlers added with PetscObjectAddOptionsHandler() */
185: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)part, PetscOptionsObject));
186: PetscOptionsEnd();
187: PetscFunctionReturn(PETSC_SUCCESS);
188: }
190: /*@
191: PetscPartitionerSetUp - Construct data structures for the PetscPartitioner
193: Collective on PetscPartitioner
195: Input Parameter:
196: . part - the PetscPartitioner object to setup
198: Level: developer
200: .seealso: `PetscPartitionerView()`, `PetscPartitionerDestroy()`
201: @*/
202: PetscErrorCode PetscPartitionerSetUp(PetscPartitioner part)
203: {
204: PetscFunctionBegin;
206: PetscTryTypeMethod(part, setup);
207: PetscFunctionReturn(PETSC_SUCCESS);
208: }
210: /*@
211: PetscPartitionerReset - Resets data structures for the PetscPartitioner
213: Collective on PetscPartitioner
215: Input Parameter:
216: . part - the PetscPartitioner object to reset
218: Level: developer
220: .seealso: `PetscPartitionerSetUp()`, `PetscPartitionerDestroy()`
221: @*/
222: PetscErrorCode PetscPartitionerReset(PetscPartitioner part)
223: {
224: PetscFunctionBegin;
226: PetscTryTypeMethod(part, reset);
227: PetscFunctionReturn(PETSC_SUCCESS);
228: }
230: /*@
231: PetscPartitionerDestroy - Destroys a PetscPartitioner object
233: Collective on PetscPartitioner
235: Input Parameter:
236: . part - the PetscPartitioner object to destroy
238: Level: developer
240: .seealso: `PetscPartitionerView()`
241: @*/
242: PetscErrorCode PetscPartitionerDestroy(PetscPartitioner *part)
243: {
244: PetscFunctionBegin;
245: if (!*part) PetscFunctionReturn(PETSC_SUCCESS);
248: if (--((PetscObject)(*part))->refct > 0) {
249: *part = NULL;
250: PetscFunctionReturn(PETSC_SUCCESS);
251: }
252: ((PetscObject)(*part))->refct = 0;
254: PetscCall(PetscPartitionerReset(*part));
256: PetscCall(PetscViewerDestroy(&(*part)->viewer));
257: PetscCall(PetscViewerDestroy(&(*part)->viewerGraph));
258: PetscTryTypeMethod((*part), destroy);
259: PetscCall(PetscHeaderDestroy(part));
260: PetscFunctionReturn(PETSC_SUCCESS);
261: }
263: /*@
264: PetscPartitionerPartition - Partition a graph
266: Collective on PetscPartitioner
268: Input Parameters:
269: + part - The PetscPartitioner
270: . nparts - Number of partitions
271: . numVertices - Number of vertices in the local part of the graph
272: . start - row pointers for the local part of the graph (CSR style)
273: . adjacency - adjacency list (CSR style)
274: . vertexSection - PetscSection describing the absolute weight of each local vertex (can be NULL)
275: - targetSection - PetscSection describing the absolute weight of each partition (can be NULL)
277: Output Parameters:
278: + partSection - The PetscSection giving the division of points by partition
279: - partition - The list of points by partition
281: Options Database:
282: + -petscpartitioner_view - View the partitioner information
283: - -petscpartitioner_view_graph - View the graph we are partitioning
285: Notes:
286: 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.
287: 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.
289: Level: developer
291: .seealso `PetscPartitionerCreate()`, `PetscPartitionerSetType()`, `PetscSectionCreate()`, `PetscSectionSetChart()`, `PetscSectionSetDof()`
292: @*/
293: PetscErrorCode PetscPartitionerPartition(PetscPartitioner part, PetscInt nparts, PetscInt numVertices, PetscInt start[], PetscInt adjacency[], PetscSection vertexSection, PetscSection targetSection, PetscSection partSection, IS *partition)
294: {
295: PetscFunctionBegin;
298: PetscCheck(nparts > 0, PetscObjectComm((PetscObject)part), PETSC_ERR_ARG_OUTOFRANGE, "Number of parts must be positive");
299: PetscCheck(numVertices >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Number of vertices must be non-negative");
300: if (numVertices && !part->noGraph) {
304: }
305: if (vertexSection) {
306: PetscInt s, e;
309: PetscCall(PetscSectionGetChart(vertexSection, &s, &e));
310: PetscCheck(s <= 0 && e >= numVertices, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Invalid vertexSection chart [%" PetscInt_FMT ",%" PetscInt_FMT ")", s, e);
311: }
312: if (targetSection) {
313: PetscInt s, e;
316: PetscCall(PetscSectionGetChart(targetSection, &s, &e));
317: PetscCheck(s <= 0 && e >= nparts, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Invalid targetSection chart [%" PetscInt_FMT ",%" PetscInt_FMT ")", s, e);
318: }
322: PetscCall(PetscSectionReset(partSection));
323: PetscCall(PetscSectionSetChart(partSection, 0, nparts));
324: if (nparts == 1) { /* quick */
325: PetscCall(PetscSectionSetDof(partSection, 0, numVertices));
326: PetscCall(ISCreateStride(PetscObjectComm((PetscObject)part), numVertices, 0, 1, partition));
327: } else PetscUseTypeMethod(part, partition, nparts, numVertices, start, adjacency, vertexSection, targetSection, partSection, partition);
328: PetscCall(PetscSectionSetUp(partSection));
329: if (part->viewerGraph) {
330: PetscViewer viewer = part->viewerGraph;
331: PetscBool isascii;
332: PetscInt v, i;
333: PetscMPIInt rank;
335: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank));
336: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
337: if (isascii) {
338: PetscCall(PetscViewerASCIIPushSynchronized(viewer));
339: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d]Nv: %" PetscInt_FMT "\n", rank, numVertices));
340: for (v = 0; v < numVertices; ++v) {
341: const PetscInt s = start[v];
342: const PetscInt e = start[v + 1];
344: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] ", rank));
345: for (i = s; i < e; ++i) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%" PetscInt_FMT " ", adjacency[i]));
346: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%" PetscInt_FMT "-%" PetscInt_FMT ")\n", s, e));
347: }
348: PetscCall(PetscViewerFlush(viewer));
349: PetscCall(PetscViewerASCIIPopSynchronized(viewer));
350: }
351: }
352: if (part->viewer) PetscCall(PetscPartitionerView(part, part->viewer));
353: PetscFunctionReturn(PETSC_SUCCESS);
354: }
356: /*@
357: PetscPartitionerCreate - Creates an empty PetscPartitioner object. The type can then be set with PetscPartitionerSetType().
359: Collective
361: Input Parameter:
362: . comm - The communicator for the PetscPartitioner object
364: Output Parameter:
365: . part - The PetscPartitioner object
367: Level: beginner
369: .seealso: `PetscPartitionerSetType()`, `PetscPartitionerDestroy()`
370: @*/
371: PetscErrorCode PetscPartitionerCreate(MPI_Comm comm, PetscPartitioner *part)
372: {
373: PetscPartitioner p;
374: const char *partitionerType = NULL;
376: PetscFunctionBegin;
378: *part = NULL;
379: PetscCall(PetscPartitionerInitializePackage());
381: PetscCall(PetscHeaderCreate(p, PETSCPARTITIONER_CLASSID, "PetscPartitioner", "Graph Partitioner", "PetscPartitioner", comm, PetscPartitionerDestroy, PetscPartitionerView));
382: PetscCall(PetscPartitionerGetDefaultType(comm, &partitionerType));
383: PetscCall(PetscPartitionerSetType(p, partitionerType));
385: p->edgeCut = 0;
386: p->balance = 0.0;
387: p->usevwgt = PETSC_TRUE;
389: *part = p;
390: PetscFunctionReturn(PETSC_SUCCESS);
391: }