Actual source code: partmatpart.c

petsc-3.14.6 2021-03-30
Report Typos and Errors
  1: #include <petscmat.h>
  2: #include <petsc/private/partitionerimpl.h>

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

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

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

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

 20:   Not Collective

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

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

 28:   Level: developer

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

212:   Level: developer

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

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

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