Actual source code: partparmetis.c
petsc-3.14.6 2021-03-30
1: #include <petsc/private/partitionerimpl.h>
3: #if defined(PETSC_HAVE_PARMETIS)
4: #include <parmetis.h>
5: #endif
7: PetscBool ParMetisPartitionerCite = PETSC_FALSE;
8: const char ParMetisPartitionerCitation[] =
9: "@article{KarypisKumar98,\n"
10: " author = {George Karypis and Vipin Kumar},\n"
11: " title = {A Parallel Algorithm for Multilevel Graph Partitioning and Sparse Matrix Ordering},\n"
12: " journal = {Journal of Parallel and Distributed Computing},\n"
13: " volume = {48},\n"
14: " pages = {71--85},\n"
15: " year = {1998}\n"
16: " doi = {https://doi.org/10.1006/jpdc.1997.1403}\n"
17: "}\n";
19: typedef struct {
20: MPI_Comm pcomm;
21: PetscInt ptype;
22: PetscReal imbalanceRatio;
23: PetscInt debugFlag;
24: PetscInt randomSeed;
25: } PetscPartitioner_ParMetis;
27: static const char *ptypes[] = {"kway", "rb"};
29: static PetscErrorCode PetscPartitionerDestroy_ParMetis(PetscPartitioner part)
30: {
31: PetscPartitioner_ParMetis *p = (PetscPartitioner_ParMetis *) part->data;
32: PetscErrorCode ierr;
35: MPI_Comm_free(&p->pcomm);
36: PetscFree(part->data);
37: return(0);
38: }
40: static PetscErrorCode PetscPartitionerView_ParMetis_ASCII(PetscPartitioner part, PetscViewer viewer)
41: {
42: PetscPartitioner_ParMetis *p = (PetscPartitioner_ParMetis *) part->data;
43: PetscErrorCode ierr;
46: PetscViewerASCIIPushTab(viewer);
47: PetscViewerASCIIPrintf(viewer, "ParMetis type: %s\n", ptypes[p->ptype]);
48: PetscViewerASCIIPrintf(viewer, "load imbalance ratio %g\n", (double) p->imbalanceRatio);
49: PetscViewerASCIIPrintf(viewer, "debug flag %D\n", p->debugFlag);
50: PetscViewerASCIIPrintf(viewer, "random seed %D\n", p->randomSeed);
51: PetscViewerASCIIPopTab(viewer);
52: return(0);
53: }
55: static PetscErrorCode PetscPartitionerView_ParMetis(PetscPartitioner part, PetscViewer viewer)
56: {
57: PetscBool iascii;
63: PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);
64: if (iascii) {PetscPartitionerView_ParMetis_ASCII(part, viewer);}
65: return(0);
66: }
68: static PetscErrorCode PetscPartitionerSetFromOptions_ParMetis(PetscOptionItems *PetscOptionsObject, PetscPartitioner part)
69: {
70: PetscPartitioner_ParMetis *p = (PetscPartitioner_ParMetis *) part->data;
71: PetscErrorCode ierr;
74: PetscOptionsHead(PetscOptionsObject, "PetscPartitioner ParMetis Options");
75: PetscOptionsEList("-petscpartitioner_parmetis_type", "Partitioning method", "", ptypes, 2, ptypes[p->ptype], &p->ptype, NULL);
76: PetscOptionsReal("-petscpartitioner_parmetis_imbalance_ratio", "Load imbalance ratio limit", "", p->imbalanceRatio, &p->imbalanceRatio, NULL);
77: PetscOptionsInt("-petscpartitioner_parmetis_debug", "Debugging flag", "", p->debugFlag, &p->debugFlag, NULL);
78: PetscOptionsInt("-petscpartitioner_parmetis_seed", "Random seed", "", p->randomSeed, &p->randomSeed, NULL);
79: PetscOptionsTail();
80: return(0);
81: }
83: static PetscErrorCode PetscPartitionerPartition_ParMetis(PetscPartitioner part, PetscInt nparts, PetscInt numVertices, PetscInt start[], PetscInt adjacency[], PetscSection vertSection, PetscSection targetSection, PetscSection partSection, IS *partition)
84: {
85: #if defined(PETSC_HAVE_PARMETIS)
86: PetscPartitioner_ParMetis *pm = (PetscPartitioner_ParMetis *) part->data;
87: MPI_Comm comm;
88: PetscInt nvtxs = numVertices; /* The number of vertices in full graph */
89: PetscInt *vtxdist; /* Distribution of vertices across processes */
90: PetscInt *xadj = start; /* Start of edge list for each vertex */
91: PetscInt *adjncy = adjacency; /* Edge lists for all vertices */
92: PetscInt *vwgt = NULL; /* Vertex weights */
93: PetscInt *adjwgt = NULL; /* Edge weights */
94: PetscInt wgtflag = 0; /* Indicates which weights are present */
95: PetscInt numflag = 0; /* Indicates initial offset (0 or 1) */
96: PetscInt ncon = 1; /* The number of weights per vertex */
97: PetscInt metis_ptype = pm->ptype; /* kway or recursive bisection */
98: real_t *tpwgts; /* The fraction of vertex weights assigned to each partition */
99: real_t *ubvec; /* The balance intolerance for vertex weights */
100: PetscInt options[64]; /* Options */
101: PetscInt v, i, *assignment, *points;
102: PetscMPIInt p, size, rank;
103: PetscBool hasempty = PETSC_FALSE;
107: PetscObjectGetComm((PetscObject) part, &comm);
108: MPI_Comm_size(comm, &size);
109: MPI_Comm_rank(comm, &rank);
110: /* Calculate vertex distribution */
111: PetscMalloc4(size+1,&vtxdist,nparts*ncon,&tpwgts,ncon,&ubvec,nvtxs,&assignment);
112: vtxdist[0] = 0;
113: MPI_Allgather(&nvtxs, 1, MPIU_INT, &vtxdist[1], 1, MPIU_INT, comm);
114: for (p = 2; p <= size; ++p) {
115: hasempty = (PetscBool)(hasempty || !vtxdist[p-1] || !vtxdist[p]);
116: vtxdist[p] += vtxdist[p-1];
117: }
118: /* null graph */
119: if (vtxdist[size] == 0) {
120: PetscFree4(vtxdist,tpwgts,ubvec,assignment);
121: ISCreateGeneral(comm, 0, NULL, PETSC_OWN_POINTER, partition);
122: return(0);
123: }
124: /* Calculate partition weights */
125: if (targetSection) {
126: PetscInt p;
127: real_t sumt = 0.0;
129: for (p = 0; p < nparts; ++p) {
130: PetscInt tpd;
132: PetscSectionGetDof(targetSection,p,&tpd);
133: sumt += tpd;
134: tpwgts[p] = tpd;
135: }
136: if (sumt) { /* METIS/ParMETIS do not like exactly zero weight */
137: for (p = 0, sumt = 0.0; p < nparts; ++p) {
138: tpwgts[p] = PetscMax(tpwgts[p],PETSC_SMALL);
139: sumt += tpwgts[p];
140: }
141: for (p = 0; p < nparts; ++p) tpwgts[p] /= sumt;
142: for (p = 0, sumt = 0.0; p < nparts-1; ++p) sumt += tpwgts[p];
143: tpwgts[nparts - 1] = 1. - sumt;
144: }
145: } else {
146: for (p = 0; p < nparts; ++p) tpwgts[p] = 1.0/nparts;
147: }
148: ubvec[0] = pm->imbalanceRatio;
150: /* Weight cells */
151: if (vertSection) {
152: PetscMalloc1(nvtxs,&vwgt);
153: for (v = 0; v < nvtxs; ++v) {
154: PetscSectionGetDof(vertSection, v, &vwgt[v]);
155: }
156: wgtflag |= 2; /* have weights on graph vertices */
157: }
159: for (p = 0; !vtxdist[p+1] && p < size; ++p);
160: if (vtxdist[p+1] == vtxdist[size]) {
161: if (rank == p) {
162: int err;
163: err = METIS_SetDefaultOptions(options); /* initialize all defaults */
164: options[METIS_OPTION_DBGLVL] = pm->debugFlag;
165: options[METIS_OPTION_SEED] = pm->randomSeed;
166: if (err != METIS_OK) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in METIS_SetDefaultOptions()");
167: if (metis_ptype == 1) {
168: PetscStackPush("METIS_PartGraphRecursive");
169: err = METIS_PartGraphRecursive(&nvtxs, &ncon, xadj, adjncy, vwgt, NULL, adjwgt, &nparts, tpwgts, ubvec, options, &part->edgeCut, assignment);
170: PetscStackPop;
171: if (err != METIS_OK) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in METIS_PartGraphRecursive()");
172: } else {
173: /*
174: It would be nice to activate the two options below, but they would need some actual testing.
175: - Turning on these options may exercise path of the METIS code that have bugs and may break production runs.
176: - If CONTIG is set to 1, METIS will exit with error if the graph is disconnected, despite the manual saying the option is ignored in such case.
177: */
178: /* options[METIS_OPTION_CONTIG] = 1; */ /* try to produce partitions that are contiguous */
179: /* options[METIS_OPTION_MINCONN] = 1; */ /* minimize the maximum degree of the subdomain graph */
180: PetscStackPush("METIS_PartGraphKway");
181: err = METIS_PartGraphKway(&nvtxs, &ncon, xadj, adjncy, vwgt, NULL, adjwgt, &nparts, tpwgts, ubvec, options, &part->edgeCut, assignment);
182: PetscStackPop;
183: if (err != METIS_OK) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in METIS_PartGraphKway()");
184: }
185: }
186: } else {
187: MPI_Comm pcomm = pm->pcomm;
189: options[0] = 1; /*use options */
190: options[1] = pm->debugFlag;
191: options[2] = (pm->randomSeed == -1) ? 15 : pm->randomSeed; /* default is GLOBAL_SEED=15 from `libparmetis/defs.h` */
193: if (hasempty) { /* parmetis does not support empty graphs on some of the processes */
194: PetscInt cnt;
196: MPI_Comm_split(pm->pcomm,!!nvtxs,rank,&pcomm);
197: for (p=0,cnt=0;p<size;p++) {
198: if (vtxdist[p+1] != vtxdist[p]) {
199: vtxdist[cnt+1] = vtxdist[p+1];
200: cnt++;
201: }
202: }
203: }
204: if (nvtxs) {
205: int err;
206: PetscStackPush("ParMETIS_V3_PartKway");
207: err = ParMETIS_V3_PartKway(vtxdist, xadj, adjncy, vwgt, adjwgt, &wgtflag, &numflag, &ncon, &nparts, tpwgts, ubvec, options, &part->edgeCut, assignment, &pcomm);
208: PetscStackPop;
209: if (err != METIS_OK) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_LIB, "Error %d in ParMETIS_V3_PartKway()", err);
210: }
211: if (hasempty) {
212: MPI_Comm_free(&pcomm);
213: }
214: }
216: /* Convert to PetscSection+IS */
217: for (v = 0; v < nvtxs; ++v) {PetscSectionAddDof(partSection, assignment[v], 1);}
218: PetscMalloc1(nvtxs, &points);
219: for (p = 0, i = 0; p < nparts; ++p) {
220: for (v = 0; v < nvtxs; ++v) {
221: if (assignment[v] == p) points[i++] = v;
222: }
223: }
224: if (i != nvtxs) SETERRQ2(comm, PETSC_ERR_PLIB, "Number of points %D should be %D", i, nvtxs);
225: ISCreateGeneral(comm, nvtxs, points, PETSC_OWN_POINTER, partition);
226: PetscFree4(vtxdist,tpwgts,ubvec,assignment);
227: PetscFree(vwgt);
228: return(0);
229: #else
230: SETERRQ(PetscObjectComm((PetscObject) part), PETSC_ERR_SUP, "Mesh partitioning needs external package support.\nPlease reconfigure with --download-parmetis.");
231: #endif
232: }
234: static PetscErrorCode PetscPartitionerInitialize_ParMetis(PetscPartitioner part)
235: {
237: part->noGraph = PETSC_FALSE;
238: part->ops->view = PetscPartitionerView_ParMetis;
239: part->ops->setfromoptions = PetscPartitionerSetFromOptions_ParMetis;
240: part->ops->destroy = PetscPartitionerDestroy_ParMetis;
241: part->ops->partition = PetscPartitionerPartition_ParMetis;
242: return(0);
243: }
245: /*MC
246: PETSCPARTITIONERPARMETIS = "parmetis" - A PetscPartitioner object using the ParMETIS library
248: Level: intermediate
250: Options Database Keys:
251: + -petscpartitioner_parmetis_type <string> - ParMETIS partitioning type. Either "kway" or "rb" (recursive bisection)
252: . -petscpartitioner_parmetis_imbalance_ratio <value> - Load imbalance ratio limit
253: . -petscpartitioner_parmetis_debug <int> - Debugging flag passed to ParMETIS/METIS routines
254: - -petscpartitioner_parmetis_seed <int> - Random seed
256: Notes: when the graph is on a single process, this partitioner actually calls METIS and not ParMETIS
258: .seealso: PetscPartitionerType, PetscPartitionerCreate(), PetscPartitionerSetType()
259: M*/
261: PETSC_EXTERN PetscErrorCode PetscPartitionerCreate_ParMetis(PetscPartitioner part)
262: {
263: PetscPartitioner_ParMetis *p;
264: PetscErrorCode ierr;
268: PetscNewLog(part, &p);
269: part->data = p;
271: MPI_Comm_dup(PetscObjectComm((PetscObject)part),&p->pcomm);
272: p->ptype = 0;
273: p->imbalanceRatio = 1.05;
274: p->debugFlag = 0;
275: p->randomSeed = -1; /* defaults to GLOBAL_SEED=15 from `libparmetis/defs.h` */
277: PetscPartitionerInitialize_ParMetis(part);
278: PetscCitationsRegister(ParMetisPartitionerCitation, &ParMetisPartitionerCite);
279: return(0);
280: }