Actual source code: partptscotch.c

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

  3: #if defined(PETSC_HAVE_PTSCOTCH)
  4: EXTERN_C_BEGIN
  5: #include <ptscotch.h>
  6: EXTERN_C_END
  7: #endif

  9: PetscBool  PTScotchPartitionerCite = PETSC_FALSE;
 10: const char PTScotchPartitionerCitation[] =
 11:   "@article{PTSCOTCH,\n"
 12:   "  author  = {C. Chevalier and F. Pellegrini},\n"
 13:   "  title   = {{PT-SCOTCH}: a tool for efficient parallel graph ordering},\n"
 14:   "  journal = {Parallel Computing},\n"
 15:   "  volume  = {34},\n"
 16:   "  number  = {6},\n"
 17:   "  pages   = {318--331},\n"
 18:   "  year    = {2008},\n"
 19:   "  doi     = {https://doi.org/10.1016/j.parco.2007.12.001}\n"
 20:   "}\n";

 22: typedef struct {
 23:   MPI_Comm  pcomm;
 24:   PetscInt  strategy;
 25:   PetscReal imbalance;
 26: } PetscPartitioner_PTScotch;

 28: #if defined(PETSC_HAVE_PTSCOTCH)

 30: #define CHKERRPTSCOTCH(ierr) do { if (ierr) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Error calling PT-Scotch library"); } while (0)

 32: static int PTScotch_Strategy(PetscInt strategy)
 33: {
 34:   switch (strategy) {
 35:   case  0: return SCOTCH_STRATDEFAULT;
 36:   case  1: return SCOTCH_STRATQUALITY;
 37:   case  2: return SCOTCH_STRATSPEED;
 38:   case  3: return SCOTCH_STRATBALANCE;
 39:   case  4: return SCOTCH_STRATSAFETY;
 40:   case  5: return SCOTCH_STRATSCALABILITY;
 41:   case  6: return SCOTCH_STRATRECURSIVE;
 42:   case  7: return SCOTCH_STRATREMAP;
 43:   default: return SCOTCH_STRATDEFAULT;
 44:   }
 45: }

 47: static PetscErrorCode PTScotch_PartGraph_Seq(SCOTCH_Num strategy, double imbalance, SCOTCH_Num n, SCOTCH_Num xadj[], SCOTCH_Num adjncy[],
 48:                                              SCOTCH_Num vtxwgt[], SCOTCH_Num adjwgt[], SCOTCH_Num nparts, SCOTCH_Num tpart[], SCOTCH_Num part[])
 49: {
 50:   SCOTCH_Graph   grafdat;
 51:   SCOTCH_Strat   stradat;
 52:   SCOTCH_Num     vertnbr = n;
 53:   SCOTCH_Num     edgenbr = xadj[n];
 54:   SCOTCH_Num*    velotab = vtxwgt;
 55:   SCOTCH_Num*    edlotab = adjwgt;
 56:   SCOTCH_Num     flagval = strategy;
 57:   double         kbalval = imbalance;

 61:   {
 62:     PetscBool flg = PETSC_TRUE;
 63:     PetscOptionsDeprecatedNoObject("-petscpartititoner_ptscotch_vertex_weight",NULL,"3.13","Use -petscpartitioner_use_vertex_weights");
 64:     PetscOptionsGetBool(NULL, NULL, "-petscpartititoner_ptscotch_vertex_weight", &flg, NULL);
 65:     if (!flg) velotab = NULL;
 66:   }
 67:   SCOTCH_graphInit(&grafdat);CHKERRPTSCOTCH(ierr);
 68:   SCOTCH_graphBuild(&grafdat, 0, vertnbr, xadj, xadj + 1, velotab, NULL, edgenbr, adjncy, edlotab);CHKERRPTSCOTCH(ierr);
 69:   SCOTCH_stratInit(&stradat);CHKERRPTSCOTCH(ierr);
 70:   SCOTCH_stratGraphMapBuild(&stradat, flagval, nparts, kbalval);CHKERRPTSCOTCH(ierr);
 71:   if (tpart) {
 72:     SCOTCH_Arch archdat;
 73:     SCOTCH_archInit(&archdat);CHKERRPTSCOTCH(ierr);
 74:     SCOTCH_archCmpltw(&archdat, nparts, tpart);CHKERRPTSCOTCH(ierr);
 75:     SCOTCH_graphMap(&grafdat, &archdat, &stradat, part);CHKERRPTSCOTCH(ierr);
 76:     SCOTCH_archExit(&archdat);
 77:   } else {
 78:     SCOTCH_graphPart(&grafdat, nparts, &stradat, part);CHKERRPTSCOTCH(ierr);
 79:   }
 80:   SCOTCH_stratExit(&stradat);
 81:   SCOTCH_graphExit(&grafdat);
 82:   return(0);
 83: }

 85: static PetscErrorCode PTScotch_PartGraph_MPI(SCOTCH_Num strategy, double imbalance, SCOTCH_Num vtxdist[], SCOTCH_Num xadj[], SCOTCH_Num adjncy[],
 86:                                              SCOTCH_Num vtxwgt[], SCOTCH_Num adjwgt[], SCOTCH_Num nparts, SCOTCH_Num tpart[], SCOTCH_Num part[], MPI_Comm comm)
 87: {
 88:   PetscMPIInt     procglbnbr;
 89:   PetscMPIInt     proclocnum;
 90:   SCOTCH_Arch     archdat;
 91:   SCOTCH_Dgraph   grafdat;
 92:   SCOTCH_Dmapping mappdat;
 93:   SCOTCH_Strat    stradat;
 94:   SCOTCH_Num      vertlocnbr;
 95:   SCOTCH_Num      edgelocnbr;
 96:   SCOTCH_Num*     veloloctab = vtxwgt;
 97:   SCOTCH_Num*     edloloctab = adjwgt;
 98:   SCOTCH_Num      flagval = strategy;
 99:   double          kbalval = imbalance;
100:   PetscErrorCode  ierr;

103:   {
104:     PetscBool flg = PETSC_TRUE;
105:     PetscOptionsDeprecatedNoObject("-petscpartititoner_ptscotch_vertex_weight",NULL,"3.13","Use -petscpartitioner_use_vertex_weights");
106:     PetscOptionsGetBool(NULL, NULL, "-petscpartititoner_ptscotch_vertex_weight", &flg, NULL);
107:     if (!flg) veloloctab = NULL;
108:   }
109:   MPI_Comm_size(comm, &procglbnbr);
110:   MPI_Comm_rank(comm, &proclocnum);
111:   vertlocnbr = vtxdist[proclocnum + 1] - vtxdist[proclocnum];
112:   edgelocnbr = xadj[vertlocnbr];

114:   SCOTCH_dgraphInit(&grafdat, comm);CHKERRPTSCOTCH(ierr);
115:   SCOTCH_dgraphBuild(&grafdat, 0, vertlocnbr, vertlocnbr, xadj, xadj + 1, veloloctab, NULL, edgelocnbr, edgelocnbr, adjncy, NULL, edloloctab);CHKERRPTSCOTCH(ierr);
116:   SCOTCH_stratInit(&stradat);CHKERRPTSCOTCH(ierr);
117:   SCOTCH_stratDgraphMapBuild(&stradat, flagval, procglbnbr, nparts, kbalval);
118:   SCOTCH_archInit(&archdat);CHKERRPTSCOTCH(ierr);
119:   if (tpart) { /* target partition weights */
120:     SCOTCH_archCmpltw(&archdat, nparts, tpart);CHKERRPTSCOTCH(ierr);
121:   } else {
122:     SCOTCH_archCmplt(&archdat, nparts);CHKERRPTSCOTCH(ierr);
123:   }
124:   SCOTCH_dgraphMapInit(&grafdat, &mappdat, &archdat, part);CHKERRPTSCOTCH(ierr);

126:   SCOTCH_dgraphMapCompute(&grafdat, &mappdat, &stradat);CHKERRPTSCOTCH(ierr);
127:   SCOTCH_dgraphMapExit(&grafdat, &mappdat);
128:   SCOTCH_archExit(&archdat);
129:   SCOTCH_stratExit(&stradat);
130:   SCOTCH_dgraphExit(&grafdat);
131:   return(0);
132: }

134: #endif /* PETSC_HAVE_PTSCOTCH */

136: static const char *const
137: PTScotchStrategyList[] = {
138:   "DEFAULT",
139:   "QUALITY",
140:   "SPEED",
141:   "BALANCE",
142:   "SAFETY",
143:   "SCALABILITY",
144:   "RECURSIVE",
145:   "REMAP"
146: };

148: static PetscErrorCode PetscPartitionerDestroy_PTScotch(PetscPartitioner part)
149: {
150:   PetscPartitioner_PTScotch *p = (PetscPartitioner_PTScotch *) part->data;
151:   PetscErrorCode             ierr;

154:   MPI_Comm_free(&p->pcomm);
155:   PetscFree(part->data);
156:   return(0);
157: }

159: static PetscErrorCode PetscPartitionerView_PTScotch_ASCII(PetscPartitioner part, PetscViewer viewer)
160: {
161:   PetscPartitioner_PTScotch *p = (PetscPartitioner_PTScotch *) part->data;
162:   PetscErrorCode            ierr;

165:   PetscViewerASCIIPushTab(viewer);
166:   PetscViewerASCIIPrintf(viewer,"using partitioning strategy %s\n",PTScotchStrategyList[p->strategy]);
167:   PetscViewerASCIIPrintf(viewer,"using load imbalance ratio %g\n",(double)p->imbalance);
168:   PetscViewerASCIIPopTab(viewer);
169:   return(0);
170: }

172: static PetscErrorCode PetscPartitionerView_PTScotch(PetscPartitioner part, PetscViewer viewer)
173: {
174:   PetscBool      iascii;

180:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);
181:   if (iascii) {PetscPartitionerView_PTScotch_ASCII(part, viewer);}
182:   return(0);
183: }

185: static PetscErrorCode PetscPartitionerSetFromOptions_PTScotch(PetscOptionItems *PetscOptionsObject, PetscPartitioner part)
186: {
187:   PetscPartitioner_PTScotch *p = (PetscPartitioner_PTScotch *) part->data;
188:   const char *const         *slist = PTScotchStrategyList;
189:   PetscInt                  nlist = (PetscInt)(sizeof(PTScotchStrategyList)/sizeof(PTScotchStrategyList[0]));
190:   PetscBool                 flag;
191:   PetscErrorCode            ierr;

194:   PetscOptionsHead(PetscOptionsObject, "PetscPartitioner PTScotch Options");
195:   PetscOptionsEList("-petscpartitioner_ptscotch_strategy","Partitioning strategy","",slist,nlist,slist[p->strategy],&p->strategy,&flag);
196:   PetscOptionsReal("-petscpartitioner_ptscotch_imbalance","Load imbalance ratio","",p->imbalance,&p->imbalance,&flag);
197:   PetscOptionsTail();
198:   return(0);
199: }

201: static PetscErrorCode PetscPartitionerPartition_PTScotch(PetscPartitioner part, PetscInt nparts, PetscInt numVertices, PetscInt start[], PetscInt adjacency[], PetscSection vertSection, PetscSection targetSection, PetscSection partSection, IS *partition)
202: {
203: #if defined(PETSC_HAVE_PTSCOTCH)
204:   MPI_Comm       comm;
205:   PetscInt       nvtxs = numVertices;   /* The number of vertices in full graph */
206:   PetscInt       *vtxdist;              /* Distribution of vertices across processes */
207:   PetscInt       *xadj   = start;       /* Start of edge list for each vertex */
208:   PetscInt       *adjncy = adjacency;   /* Edge lists for all vertices */
209:   PetscInt       *vwgt   = NULL;        /* Vertex weights */
210:   PetscInt       *adjwgt = NULL;        /* Edge weights */
211:   PetscInt       v, i, *assignment, *points;
212:   PetscMPIInt    size, rank, p;
213:   PetscBool      hasempty = PETSC_FALSE;
214:   PetscInt       *tpwgts = NULL;

218:   PetscObjectGetComm((PetscObject)part,&comm);
219:   MPI_Comm_size(comm, &size);
220:   MPI_Comm_rank(comm, &rank);
221:   PetscMalloc2(size+1,&vtxdist,PetscMax(nvtxs,1),&assignment);
222:   /* Calculate vertex distribution */
223:   vtxdist[0] = 0;
224:   MPI_Allgather(&nvtxs, 1, MPIU_INT, &vtxdist[1], 1, MPIU_INT, comm);
225:   for (p = 2; p <= size; ++p) {
226:     hasempty = (PetscBool)(hasempty || !vtxdist[p-1] || !vtxdist[p]);
227:     vtxdist[p] += vtxdist[p-1];
228:   }
229:   /* null graph */
230:   if (vtxdist[size] == 0) {
231:     PetscFree2(vtxdist, assignment);
232:     ISCreateGeneral(comm, 0, NULL, PETSC_OWN_POINTER, partition);
233:     return(0);
234:   }

236:   /* Calculate vertex weights */
237:   if (vertSection) {
238:     PetscMalloc1(nvtxs,&vwgt);
239:     for (v = 0; v < nvtxs; ++v) {
240:       PetscSectionGetDof(vertSection, v, &vwgt[v]);
241:     }
242:   }

244:   /* Calculate partition weights */
245:   if (targetSection) {
246:     PetscInt sumw;

248:     PetscCalloc1(nparts,&tpwgts);
249:     for (p = 0, sumw = 0; p < nparts; ++p) {
250:       PetscSectionGetDof(targetSection,p,&tpwgts[p]);
251:       sumw += tpwgts[p];
252:     }
253:     if (!sumw) {
254:       PetscFree(tpwgts);
255:     }
256:   }

258:   {
259:     PetscPartitioner_PTScotch *pts = (PetscPartitioner_PTScotch *) part->data;
260:     int                       strat = PTScotch_Strategy(pts->strategy);
261:     double                    imbal = (double)pts->imbalance;

263:     for (p = 0; !vtxdist[p+1] && p < size; ++p);
264:     if (vtxdist[p+1] == vtxdist[size]) {
265:       if (rank == p) {
266:         PTScotch_PartGraph_Seq(strat, imbal, nvtxs, xadj, adjncy, vwgt, adjwgt, nparts, tpwgts, assignment);
267:       }
268:     } else {
269:       MPI_Comm pcomm = pts->pcomm;

271:       if (hasempty) {
272:         PetscInt cnt;

274:         MPI_Comm_split(pts->pcomm,!!nvtxs,rank,&pcomm);
275:         for (p=0,cnt=0;p<size;p++) {
276:           if (vtxdist[p+1] != vtxdist[p]) {
277:             vtxdist[cnt+1] = vtxdist[p+1];
278:             cnt++;
279:           }
280:         }
281:       };
282:       if (nvtxs) {
283:         PTScotch_PartGraph_MPI(strat, imbal, vtxdist, xadj, adjncy, vwgt, adjwgt, nparts, tpwgts, assignment, pcomm);
284:       }
285:       if (hasempty) {
286:         MPI_Comm_free(&pcomm);
287:       }
288:     }
289:   }
290:   PetscFree(vwgt);
291:   PetscFree(tpwgts);

293:   /* Convert to PetscSection+IS */
294:   for (v = 0; v < nvtxs; ++v) {PetscSectionAddDof(partSection, assignment[v], 1);}
295:   PetscMalloc1(nvtxs, &points);
296:   for (p = 0, i = 0; p < nparts; ++p) {
297:     for (v = 0; v < nvtxs; ++v) {
298:       if (assignment[v] == p) points[i++] = v;
299:     }
300:   }
301:   if (i != nvtxs) SETERRQ2(comm, PETSC_ERR_PLIB, "Number of points %D should be %D", i, nvtxs);
302:   ISCreateGeneral(comm, nvtxs, points, PETSC_OWN_POINTER, partition);

304:   PetscFree2(vtxdist,assignment);
305:   return(0);
306: #else
307:   SETERRQ(PetscObjectComm((PetscObject) part), PETSC_ERR_SUP, "Mesh partitioning needs external package support.\nPlease reconfigure with --download-ptscotch.");
308: #endif
309: }

311: static PetscErrorCode PetscPartitionerInitialize_PTScotch(PetscPartitioner part)
312: {
314:   part->noGraph             = PETSC_FALSE;
315:   part->ops->view           = PetscPartitionerView_PTScotch;
316:   part->ops->destroy        = PetscPartitionerDestroy_PTScotch;
317:   part->ops->partition      = PetscPartitionerPartition_PTScotch;
318:   part->ops->setfromoptions = PetscPartitionerSetFromOptions_PTScotch;
319:   return(0);
320: }

322: /*MC
323:   PETSCPARTITIONERPTSCOTCH = "ptscotch" - A PetscPartitioner object using the PT-Scotch library

325:   Level: intermediate

327:   Options Database Keys:
328: +  -petscpartitioner_ptscotch_strategy <string> - PT-Scotch strategy. Choose one of default quality speed balance safety scalability recursive remap
329: -  -petscpartitioner_ptscotch_imbalance <val> - Load imbalance ratio

331:   Notes: when the graph is on a single process, this partitioner actually uses Scotch and not PT-Scotch

333: .seealso: PetscPartitionerType, PetscPartitionerCreate(), PetscPartitionerSetType()
334: M*/

336: PETSC_EXTERN PetscErrorCode PetscPartitionerCreate_PTScotch(PetscPartitioner part)
337: {
338:   PetscPartitioner_PTScotch *p;
339:   PetscErrorCode          ierr;

343:   PetscNewLog(part, &p);
344:   part->data = p;

346:   MPI_Comm_dup(PetscObjectComm((PetscObject)part),&p->pcomm);
347:   p->strategy  = 0;
348:   p->imbalance = 0.01;

350:   PetscPartitionerInitialize_PTScotch(part);
351:   PetscCitationsRegister(PTScotchPartitionerCitation, &PTScotchPartitionerCite);
352:   return(0);
353: }