Actual source code: partsimple.c

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

  4: typedef struct {
  5:   PetscInt dummy;
  6: } PetscPartitioner_Simple;

  8: static PetscErrorCode PetscPartitionerDestroy_Simple(PetscPartitioner part)
  9: {

 13:   PetscFree(part->data);
 14:   return(0);
 15: }

 17: static PetscErrorCode PetscPartitionerView_Simple_ASCII(PetscPartitioner part, PetscViewer viewer)
 18: {
 20:   return(0);
 21: }

 23: static PetscErrorCode PetscPartitionerView_Simple(PetscPartitioner part, PetscViewer viewer)
 24: {
 25:   PetscBool      iascii;

 31:   PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);
 32:   if (iascii) {PetscPartitionerView_Simple_ASCII(part, viewer);}
 33:   return(0);
 34: }

 36: static PetscErrorCode PetscPartitionerPartition_Simple(PetscPartitioner part, PetscInt nparts, PetscInt numVertices, PetscInt start[], PetscInt adjacency[], PetscSection vertSection, PetscSection targetSection, PetscSection partSection, IS *partition)
 37: {
 38:   MPI_Comm       comm;
 39:   PetscInt       np, *tpwgts = NULL, sumw = 0, numVerticesGlobal  = 0;
 40:   PetscMPIInt    size;

 44:   if (vertSection) { PetscInfo(part,"PETSCPARTITIONERSIMPLE ignores vertex weights\n"); }
 45:   comm = PetscObjectComm((PetscObject)part);
 46:   MPI_Comm_size(comm,&size);
 47:   if (targetSection) {
 48:     MPIU_Allreduce(&numVertices, &numVerticesGlobal, 1, MPIU_INT, MPI_SUM, comm);
 49:     PetscCalloc1(nparts,&tpwgts);
 50:     for (np = 0; np < nparts; ++np) {
 51:       PetscSectionGetDof(targetSection,np,&tpwgts[np]);
 52:       sumw += tpwgts[np];
 53:     }
 54:     if (!sumw) {
 55:       PetscFree(tpwgts);
 56:     } else {
 57:       PetscInt m,mp;
 58:       for (np = 0; np < nparts; ++np) tpwgts[np] = (tpwgts[np]*numVerticesGlobal)/sumw;
 59:       for (np = 0, m = -1, mp = 0, sumw = 0; np < nparts; ++np) {
 60:         if (m < tpwgts[np]) { m = tpwgts[np]; mp = np; }
 61:         sumw += tpwgts[np];
 62:       }
 63:       if (sumw != numVerticesGlobal) tpwgts[mp] += numVerticesGlobal - sumw;
 64:     }
 65:   }

 67:   ISCreateStride(PETSC_COMM_SELF, numVertices, 0, 1, partition);
 68:   if (size == 1) {
 69:     if (tpwgts) {
 70:       for (np = 0; np < nparts; ++np) {
 71:         PetscSectionSetDof(partSection, np, tpwgts[np]);
 72:       }
 73:     } else {
 74:       for (np = 0; np < nparts; ++np) {
 75:         PetscSectionSetDof(partSection, np, numVertices/nparts + ((numVertices % nparts) > np));
 76:       }
 77:     }
 78:   } else {
 79:     if (tpwgts) {
 80:       Vec         v;
 81:       PetscScalar *array;
 82:       PetscInt    st,j;
 83:       PetscMPIInt rank;

 85:       VecCreate(comm,&v);
 86:       VecSetSizes(v,numVertices,numVerticesGlobal);
 87:       VecSetType(v,VECSTANDARD);
 88:       MPI_Comm_rank(comm,&rank);
 89:       for (np = 0,st = 0; np < nparts; ++np) {
 90:         if (rank == np || (rank == size-1 && size < nparts && np >= size)) {
 91:           for (j = 0; j < tpwgts[np]; j++) {
 92:             VecSetValue(v,st+j,np,INSERT_VALUES);
 93:           }
 94:         }
 95:         st += tpwgts[np];
 96:       }
 97:       VecAssemblyBegin(v);
 98:       VecAssemblyEnd(v);
 99:       VecGetArray(v,&array);
100:       for (j = 0; j < numVertices; ++j) {
101:         PetscSectionAddDof(partSection,PetscRealPart(array[j]),1);
102:       }
103:       VecRestoreArray(v,&array);
104:       VecDestroy(&v);
105:     } else {
106:       PetscMPIInt rank;
107:       PetscInt nvGlobal, *offsets, myFirst, myLast;

109:       PetscMalloc1(size+1,&offsets);
110:       offsets[0] = 0;
111:       MPI_Allgather(&numVertices,1,MPIU_INT,&offsets[1],1,MPIU_INT,comm);
112:       for (np = 2; np <= size; np++) {
113:         offsets[np] += offsets[np-1];
114:       }
115:       nvGlobal = offsets[size];
116:       MPI_Comm_rank(comm,&rank);
117:       myFirst = offsets[rank];
118:       myLast  = offsets[rank + 1] - 1;
119:       PetscFree(offsets);
120:       if (numVertices) {
121:         PetscInt firstPart = 0, firstLargePart = 0;
122:         PetscInt lastPart = 0, lastLargePart = 0;
123:         PetscInt rem = nvGlobal % nparts;
124:         PetscInt pSmall = nvGlobal/nparts;
125:         PetscInt pBig = nvGlobal/nparts + 1;

127:         if (rem) {
128:           firstLargePart = myFirst / pBig;
129:           lastLargePart  = myLast  / pBig;

131:           if (firstLargePart < rem) {
132:             firstPart = firstLargePart;
133:           } else {
134:             firstPart = rem + (myFirst - (rem * pBig)) / pSmall;
135:           }
136:           if (lastLargePart < rem) {
137:             lastPart = lastLargePart;
138:           } else {
139:             lastPart = rem + (myLast - (rem * pBig)) / pSmall;
140:           }
141:         } else {
142:           firstPart = myFirst / (nvGlobal/nparts);
143:           lastPart  = myLast  / (nvGlobal/nparts);
144:         }

146:         for (np = firstPart; np <= lastPart; np++) {
147:           PetscInt PartStart =  np    * (nvGlobal/nparts) + PetscMin(nvGlobal % nparts,np);
148:           PetscInt PartEnd   = (np+1) * (nvGlobal/nparts) + PetscMin(nvGlobal % nparts,np+1);

150:           PartStart = PetscMax(PartStart,myFirst);
151:           PartEnd   = PetscMin(PartEnd,myLast+1);
152:           PetscSectionSetDof(partSection,np,PartEnd-PartStart);
153:         }
154:       }
155:     }
156:   }
157:   PetscFree(tpwgts);
158:   return(0);
159: }

161: static PetscErrorCode PetscPartitionerInitialize_Simple(PetscPartitioner part)
162: {
164:   part->noGraph        = PETSC_TRUE;
165:   part->ops->view      = PetscPartitionerView_Simple;
166:   part->ops->destroy   = PetscPartitionerDestroy_Simple;
167:   part->ops->partition = PetscPartitionerPartition_Simple;
168:   return(0);
169: }

171: /*MC
172:   PETSCPARTITIONERSIMPLE = "simple" - A PetscPartitioner object

174:   Level: intermediate

176: .seealso: PetscPartitionerType, PetscPartitionerCreate(), PetscPartitionerSetType()
177: M*/

179: PETSC_EXTERN PetscErrorCode PetscPartitionerCreate_Simple(PetscPartitioner part)
180: {
181:   PetscPartitioner_Simple *p;
182:   PetscErrorCode           ierr;

186:   PetscNewLog(part, &p);
187:   part->data = p;

189:   PetscPartitionerInitialize_Simple(part);
190:   return(0);
191: }