Actual source code: partsimple.c
petsc-master 2020-08-25
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: }