Actual source code: petscpartmatpart.c

petsc-master 2020-05-29
Report Typos and Errors
  1:  #include <petsc/private/dmpleximpl.h>

  3: typedef struct {
  4:   MatPartitioning mp;
  5: } PetscPartitioner_MatPartitioning;

  7: static PetscErrorCode PetscPartitionerMatPartitioningGetMatPartitioning_MatPartitioning(PetscPartitioner part, MatPartitioning *mp)
  8: {
  9:   PetscPartitioner_MatPartitioning  *p = (PetscPartitioner_MatPartitioning *) part->data;

 12:   *mp = p->mp;
 13:   return(0);
 14: }

 16: /*@C
 17:   PetscPartitionerMatPartitioningGetMatPartitioning - Get a MatPartitioning instance wrapped by this PetscPartitioner.

 19:   Not Collective

 21:   Input Parameters:
 22: . part     - The PetscPartitioner

 24:   Output Parameters:
 25: . mp       - The MatPartitioning

 27:   Level: developer

 29: .seealso PetscPartitionerMatPartitioningSetMatPartitioning(), DMPlexDistribute(), PetscPartitionerCreate()
 30: @*/
 31: PetscErrorCode PetscPartitionerMatPartitioningGetMatPartitioning(PetscPartitioner part, MatPartitioning *mp)
 32: {
 33:   PetscErrorCode          ierr;

 38:   PetscUseMethod(part,"PetscPartitionerMatPartitioningGetMatPartitioning_C",(PetscPartitioner,MatPartitioning*),(part,mp));
 39:   return(0);
 40: }

 42: static PetscErrorCode PetscPartitionerDestroy_MatPartitioning(PetscPartitioner part)
 43: {
 44:   PetscPartitioner_MatPartitioning  *p = (PetscPartitioner_MatPartitioning *) part->data;
 45:   PetscErrorCode                    ierr;

 48:   MatPartitioningDestroy(&p->mp);
 49:   PetscFree(p);
 50:   return(0);
 51: }

 53: static PetscErrorCode PetscPartitionerView_MatPartitioning_Ascii(PetscPartitioner part, PetscViewer viewer)
 54: {
 55:   PetscPartitioner_MatPartitioning  *p = (PetscPartitioner_MatPartitioning *) part->data;
 56:   PetscViewerFormat                 format;
 57:   PetscErrorCode                    ierr;

 60:   PetscViewerGetFormat(viewer, &format);
 61:   PetscViewerASCIIPrintf(viewer, "MatPartitioning Graph Partitioner:\n");
 62:   PetscViewerASCIIPushTab(viewer);
 63:   if (p->mp) {MatPartitioningView(p->mp, viewer);}
 64:   PetscViewerASCIIPopTab(viewer);
 65:   return(0);
 66: }

 68: static PetscErrorCode PetscPartitionerView_MatPartitioning(PetscPartitioner part, PetscViewer viewer)
 69: {
 70:   PetscBool      iascii;

 76:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);
 77:   if (iascii) {PetscPartitionerView_MatPartitioning_Ascii(part, viewer);}
 78:   return(0);
 79: }

 81: static PetscErrorCode PetscPartitionerSetFromOptions_MatPartitioning(PetscOptionItems *PetscOptionsObject, PetscPartitioner part)
 82: {
 83:   PetscPartitioner_MatPartitioning  *p = (PetscPartitioner_MatPartitioning *) part->data;
 84:   PetscErrorCode                    ierr;

 87:   PetscObjectSetOptionsPrefix((PetscObject)p->mp,((PetscObject)part)->prefix);
 88:   MatPartitioningSetFromOptions(p->mp);
 89:   return(0);
 90: }

 92: static PetscErrorCode PetscPartitionerPartition_MatPartitioning(PetscPartitioner part, PetscInt nparts, PetscInt numVertices, PetscInt start[], PetscInt adjacency[], PetscSection vertSection, PetscSection targetSection, PetscSection partSection, IS *is)
 93: {
 94:   PetscPartitioner_MatPartitioning  *p = (PetscPartitioner_MatPartitioning *) part->data;
 95:   Mat                               matadj;
 96:   IS                                is1, is2, is3;
 97:   PetscReal                         *tpwgts = NULL;
 98:   PetscInt                          numVerticesGlobal, numEdges;
 99:   PetscInt                          *i, *j, *vwgt = NULL;
100:   MPI_Comm                          comm;
101:   PetscErrorCode                    ierr;

104:   PetscObjectGetComm((PetscObject)part, &comm);

106:   /* TODO: MatCreateMPIAdj should maybe take global number of ROWS */
107:   /* TODO: And vertex distribution in PetscPartitionerPartition_ParMetis should be done using PetscSplitOwnership */
108:   numVerticesGlobal = PETSC_DECIDE;
109:   PetscSplitOwnership(comm, &numVertices, &numVerticesGlobal);

111:   /* copy arrays to avoid memory errors because MatMPIAdjSetPreallocation copies just pointers */
112:   numEdges = start[numVertices];
113:   PetscMalloc1(numVertices+1, &i);
114:   PetscMalloc1(numEdges, &j);
115:   PetscArraycpy(i, start, numVertices+1);
116:   PetscArraycpy(j, adjacency, numEdges);

118:   /* construct the adjacency matrix */
119:   MatCreateMPIAdj(comm, numVertices, numVerticesGlobal, i, j, NULL, &matadj);
120:   MatPartitioningSetAdjacency(p->mp, matadj);
121:   MatPartitioningSetNParts(p->mp, nparts);

123:   /* calculate partition weights */
124:   if (targetSection) {
125:     PetscReal sumt;
126:     PetscInt  p;

128:     sumt = 0.0;
129:     PetscMalloc1(nparts,&tpwgts);
130:     for (p = 0; p < nparts; ++p) {
131:       PetscInt tpd;

133:       PetscSectionGetDof(targetSection,p,&tpd);
134:       sumt += tpd;
135:       tpwgts[p] = tpd;
136:     }
137:     if (sumt) { /* METIS/ParMETIS do not like exactly zero weight */
138:       for (p = 0, sumt = 0.0; p < nparts; ++p) {
139:         tpwgts[p] = PetscMax(tpwgts[p],PETSC_SMALL);
140:         sumt += tpwgts[p];
141:       }
142:       for (p = 0; p < nparts; ++p) tpwgts[p] /= sumt;
143:       for (p = 0, sumt = 0.0; p < nparts-1; ++p) sumt += tpwgts[p];
144:       tpwgts[nparts - 1] = 1. - sumt;
145:     } else {
146:       PetscFree(tpwgts);
147:     }
148:   }
149:   MatPartitioningSetPartitionWeights(p->mp, tpwgts);

151:   /* calculate vertex weights */
152:   if (vertSection) {
153:     PetscInt v;

155:     PetscMalloc1(numVertices,&vwgt);
156:     for (v = 0; v < numVertices; ++v) {
157:       PetscSectionGetDof(vertSection, v, &vwgt[v]);
158:     }
159:   }
160:   MatPartitioningSetVertexWeights(p->mp, vwgt);

162:   /* apply the partitioning */
163:   MatPartitioningApply(p->mp, &is1);

165:   /* construct the PetscSection */
166:   {
167:     PetscInt v;
168:     const PetscInt *assignment_arr;

170:     ISGetIndices(is1, &assignment_arr);
171:     for (v = 0; v < numVertices; ++v) {PetscSectionAddDof(partSection, assignment_arr[v], 1);}
172:     ISRestoreIndices(is1, &assignment_arr);
173:   }

175:   /* convert assignment IS to global numbering IS */
176:   ISPartitioningToNumbering(is1, &is2);
177:   ISDestroy(&is1);

179:   /* renumber IS into local numbering */
180:   ISOnComm(is2, PETSC_COMM_SELF, PETSC_USE_POINTER, &is1);
181:   ISRenumber(is1, NULL, NULL, &is3);
182:   ISDestroy(&is1);
183:   ISDestroy(&is2);

185:   /* invert IS */
186:   ISSetPermutation(is3);
187:   ISInvertPermutation(is3, numVertices, &is1);
188:   ISDestroy(&is3);

190:   MatDestroy(&matadj);
191:   *is = is1;
192:   return(0);
193: }

195: static PetscErrorCode PetscPartitionerInitialize_MatPartitioning(PetscPartitioner part)
196: {
197:   PetscErrorCode                    ierr;

200:   part->ops->view           = PetscPartitionerView_MatPartitioning;
201:   part->ops->setfromoptions = PetscPartitionerSetFromOptions_MatPartitioning;
202:   part->ops->destroy        = PetscPartitionerDestroy_MatPartitioning;
203:   part->ops->partition      = PetscPartitionerPartition_MatPartitioning;
204:   PetscObjectComposeFunction((PetscObject)part,"PetscPartitionerMatPartitioningGetMatPartitioning_C",PetscPartitionerMatPartitioningGetMatPartitioning_MatPartitioning);
205:   return(0);
206: }

208: /*MC
209:   PETSCPARTITIONERMATPARTITIONING = "matpartitioning" - A PetscPartitioner object

211:   Level: developer

213: .seealso: PetscPartitionerType, PetscPartitionerCreate(), PetscPartitionerSetType()
214: M*/

216: PETSC_EXTERN PetscErrorCode PetscPartitionerCreate_MatPartitioning(PetscPartitioner part)
217: {
218:   PetscPartitioner_MatPartitioning  *p;
219:   PetscErrorCode                    ierr;

223:   PetscNewLog(part, &p);
224:   part->data = p;
225:   PetscPartitionerInitialize_MatPartitioning(part);
226:   MatPartitioningCreate(PetscObjectComm((PetscObject)part), &p->mp);
227:   return(0);
228: }