Actual source code: partmatpart.c

  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;

 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 `DMPlexDistribute()`, `PetscPartitionerCreate()`
 30: @*/
 31: PetscErrorCode PetscPartitionerMatPartitioningGetMatPartitioning(PetscPartitioner part, MatPartitioning *mp)
 32: {
 35:   PetscUseMethod(part, "PetscPartitionerMatPartitioningGetMatPartitioning_C", (PetscPartitioner, MatPartitioning *), (part, mp));
 36:   return 0;
 37: }

 39: static PetscErrorCode PetscPartitionerDestroy_MatPartitioning(PetscPartitioner part)
 40: {
 41:   PetscPartitioner_MatPartitioning *p = (PetscPartitioner_MatPartitioning *)part->data;

 43:   MatPartitioningDestroy(&p->mp);
 44:   PetscObjectComposeFunction((PetscObject)part, "PetscPartitionerMatPartitioningGetMatPartitioning_C", NULL);
 45:   PetscFree(part->data);
 46:   return 0;
 47: }

 49: static PetscErrorCode PetscPartitionerView_MatPartitioning_ASCII(PetscPartitioner part, PetscViewer viewer)
 50: {
 51:   PetscPartitioner_MatPartitioning *p = (PetscPartitioner_MatPartitioning *)part->data;
 52:   PetscViewerFormat                 format;

 54:   PetscViewerGetFormat(viewer, &format);
 55:   PetscViewerASCIIPrintf(viewer, "MatPartitioning Graph Partitioner:\n");
 56:   PetscViewerASCIIPushTab(viewer);
 57:   if (p->mp) MatPartitioningView(p->mp, viewer);
 58:   PetscViewerASCIIPopTab(viewer);
 59:   return 0;
 60: }

 62: static PetscErrorCode PetscPartitionerView_MatPartitioning(PetscPartitioner part, PetscViewer viewer)
 63: {
 64:   PetscBool iascii;

 68:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
 69:   if (iascii) PetscPartitionerView_MatPartitioning_ASCII(part, viewer);
 70:   return 0;
 71: }

 73: static PetscErrorCode PetscPartitionerSetFromOptions_MatPartitioning(PetscPartitioner part, PetscOptionItems *PetscOptionsObject)
 74: {
 75:   PetscPartitioner_MatPartitioning *p = (PetscPartitioner_MatPartitioning *)part->data;

 77:   PetscObjectSetOptionsPrefix((PetscObject)p->mp, ((PetscObject)part)->prefix);
 78:   MatPartitioningSetFromOptions(p->mp);
 79:   return 0;
 80: }

 82: static PetscErrorCode PetscPartitionerPartition_MatPartitioning(PetscPartitioner part, PetscInt nparts, PetscInt numVertices, PetscInt start[], PetscInt adjacency[], PetscSection vertSection, PetscSection targetSection, PetscSection partSection, IS *is)
 83: {
 84:   PetscPartitioner_MatPartitioning *p = (PetscPartitioner_MatPartitioning *)part->data;
 85:   Mat                               matadj;
 86:   IS                                is1, is2, is3;
 87:   PetscReal                        *tpwgts = NULL;
 88:   PetscInt                          numVerticesGlobal, numEdges;
 89:   PetscInt                         *i, *j, *vwgt = NULL;
 90:   MPI_Comm                          comm;

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

 94:   /* TODO: MatCreateMPIAdj should maybe take global number of ROWS */
 95:   /* TODO: And vertex distribution in PetscPartitionerPartition_ParMetis should be done using PetscSplitOwnership */
 96:   numVerticesGlobal = PETSC_DECIDE;
 97:   PetscSplitOwnership(comm, &numVertices, &numVerticesGlobal);

 99:   /* copy arrays to avoid memory errors because MatMPIAdjSetPreallocation copies just pointers */
100:   numEdges = start[numVertices];
101:   PetscMalloc1(numVertices + 1, &i);
102:   PetscMalloc1(numEdges, &j);
103:   PetscArraycpy(i, start, numVertices + 1);
104:   PetscArraycpy(j, adjacency, numEdges);

106:   /* construct the adjacency matrix */
107:   MatCreateMPIAdj(comm, numVertices, numVerticesGlobal, i, j, NULL, &matadj);
108:   MatPartitioningSetAdjacency(p->mp, matadj);
109:   MatPartitioningSetNParts(p->mp, nparts);

111:   /* calculate partition weights */
112:   if (targetSection) {
113:     PetscReal sumt;
114:     PetscInt  p;

116:     sumt = 0.0;
117:     PetscMalloc1(nparts, &tpwgts);
118:     for (p = 0; p < nparts; ++p) {
119:       PetscInt tpd;

121:       PetscSectionGetDof(targetSection, p, &tpd);
122:       sumt += tpd;
123:       tpwgts[p] = tpd;
124:     }
125:     if (sumt) { /* METIS/ParMETIS do not like exactly zero weight */
126:       for (p = 0, sumt = 0.0; p < nparts; ++p) {
127:         tpwgts[p] = PetscMax(tpwgts[p], PETSC_SMALL);
128:         sumt += tpwgts[p];
129:       }
130:       for (p = 0; p < nparts; ++p) tpwgts[p] /= sumt;
131:       for (p = 0, sumt = 0.0; p < nparts - 1; ++p) sumt += tpwgts[p];
132:       tpwgts[nparts - 1] = 1. - sumt;
133:     } else {
134:       PetscFree(tpwgts);
135:     }
136:   }
137:   MatPartitioningSetPartitionWeights(p->mp, tpwgts);

139:   /* calculate vertex weights */
140:   if (vertSection) {
141:     PetscInt v;

143:     PetscMalloc1(numVertices, &vwgt);
144:     for (v = 0; v < numVertices; ++v) PetscSectionGetDof(vertSection, v, &vwgt[v]);
145:   }
146:   MatPartitioningSetVertexWeights(p->mp, vwgt);

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

151:   /* construct the PetscSection */
152:   {
153:     PetscInt        v;
154:     const PetscInt *assignment_arr;

156:     ISGetIndices(is1, &assignment_arr);
157:     for (v = 0; v < numVertices; ++v) PetscSectionAddDof(partSection, assignment_arr[v], 1);
158:     ISRestoreIndices(is1, &assignment_arr);
159:   }

161:   /* convert assignment IS to global numbering IS */
162:   ISPartitioningToNumbering(is1, &is2);
163:   ISDestroy(&is1);

165:   /* renumber IS into local numbering */
166:   ISOnComm(is2, PETSC_COMM_SELF, PETSC_USE_POINTER, &is1);
167:   ISRenumber(is1, NULL, NULL, &is3);
168:   ISDestroy(&is1);
169:   ISDestroy(&is2);

171:   /* invert IS */
172:   ISSetPermutation(is3);
173:   ISInvertPermutation(is3, numVertices, &is1);
174:   ISDestroy(&is3);

176:   MatDestroy(&matadj);
177:   *is = is1;
178:   return 0;
179: }

181: static PetscErrorCode PetscPartitionerInitialize_MatPartitioning(PetscPartitioner part)
182: {
183:   part->ops->view           = PetscPartitionerView_MatPartitioning;
184:   part->ops->setfromoptions = PetscPartitionerSetFromOptions_MatPartitioning;
185:   part->ops->destroy        = PetscPartitionerDestroy_MatPartitioning;
186:   part->ops->partition      = PetscPartitionerPartition_MatPartitioning;
187:   PetscObjectComposeFunction((PetscObject)part, "PetscPartitionerMatPartitioningGetMatPartitioning_C", PetscPartitionerMatPartitioningGetMatPartitioning_MatPartitioning);
188:   return 0;
189: }

191: /*MC
192:   PETSCPARTITIONERMATPARTITIONING = "matpartitioning" - A PetscPartitioner object

194:   Level: developer

196: .seealso: `PetscPartitionerType`, `PetscPartitionerCreate()`, `PetscPartitionerSetType()`
197: M*/

199: PETSC_EXTERN PetscErrorCode PetscPartitionerCreate_MatPartitioning(PetscPartitioner part)
200: {
201:   PetscPartitioner_MatPartitioning *p;

204:   PetscNew(&p);
205:   part->data = p;
206:   PetscPartitionerInitialize_MatPartitioning(part);
207:   MatPartitioningCreate(PetscObjectComm((PetscObject)part), &p->mp);
208:   return 0;
209: }