Actual source code: ex18.c

  1: static char help[] = "Test PetscSFConcatenate()\n\n";

  3: #include <petscsf.h>

  5: typedef struct {
  6:   MPI_Comm                   comm;
  7:   PetscMPIInt                rank, size;
  8:   PetscInt                   leaveStep, nsfs, n;
  9:   PetscBool                  sparseLeaves;
 10:   PetscBool                  compare;
 11:   PetscBool                  irregular;
 12:   PetscSFConcatenateRootMode rootMode;
 13:   PetscViewer                viewer;
 14: } AppCtx;

 16: static PetscErrorCode GetOptions(MPI_Comm comm, AppCtx *ctx)
 17: {
 18:   PetscViewerFormat format;

 20:   PetscFunctionBegin;
 21:   ctx->comm         = comm;
 22:   ctx->nsfs         = 3;
 23:   ctx->n            = 1;
 24:   ctx->leaveStep    = 1;
 25:   ctx->sparseLeaves = PETSC_FALSE;
 26:   ctx->compare      = PETSC_FALSE;
 27:   ctx->irregular    = PETSC_FALSE;
 28:   ctx->rootMode     = PETSCSF_CONCATENATE_ROOTMODE_LOCAL;
 29:   ctx->viewer       = NULL;
 30:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-nsfs", &ctx->nsfs, NULL));
 31:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &ctx->n, NULL));
 32:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-leave_step", &ctx->leaveStep, NULL));
 33:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-irregular", &ctx->irregular, NULL));
 34:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-compare_to_reference", &ctx->compare, NULL));
 35:   PetscCall(PetscOptionsGetEnum(NULL, NULL, "-root_mode", PetscSFConcatenateRootModes, (PetscEnum *)&ctx->rootMode, NULL));
 36:   PetscCall(PetscOptionsCreateViewer(comm, NULL, NULL, "-sf_view", &ctx->viewer, &format, NULL));
 37:   if (ctx->viewer) PetscCall(PetscViewerPushFormat(ctx->viewer, format));
 38:   ctx->sparseLeaves = (PetscBool)(ctx->leaveStep != 1);
 39:   PetscCallMPI(MPI_Comm_size(comm, &ctx->size));
 40:   PetscCallMPI(MPI_Comm_rank(comm, &ctx->rank));
 41:   PetscFunctionReturn(PETSC_SUCCESS);
 42: }

 44: static PetscErrorCode PetscSFCheckEqual_Private(PetscSF sf0, PetscSF sf1)
 45: {
 46:   PetscInt  nRoot, nLeave;
 47:   Vec       vecRoot0, vecLeave0, vecRoot1, vecLeave1;
 48:   MPI_Comm  comm;
 49:   PetscBool flg;

 51:   PetscFunctionBegin;
 52:   PetscCall(PetscObjectGetComm((PetscObject)sf0, &comm));
 53:   PetscCall(PetscSFGetGraph(sf0, &nRoot, NULL, NULL, NULL));
 54:   PetscCall(PetscSFGetLeafRange(sf0, NULL, &nLeave));
 55:   nLeave++;
 56:   PetscCall(VecCreateFromOptions(comm, NULL, 1, nRoot, PETSC_DECIDE, &vecRoot0));
 57:   PetscCall(VecCreateFromOptions(comm, NULL, 1, nLeave, PETSC_DECIDE, &vecLeave0));
 58:   PetscCall(VecDuplicate(vecRoot0, &vecRoot1));
 59:   PetscCall(VecDuplicate(vecLeave0, &vecLeave1));
 60:   {
 61:     PetscRandom rand;

 63:     PetscCall(PetscRandomCreate(comm, &rand));
 64:     PetscCall(PetscRandomSetFromOptions(rand));
 65:     PetscCall(VecSetRandom(vecRoot0, rand));
 66:     PetscCall(VecSetRandom(vecLeave0, rand));
 67:     PetscCall(VecCopy(vecRoot0, vecRoot1));
 68:     PetscCall(VecCopy(vecLeave0, vecLeave1));
 69:     PetscCall(PetscRandomDestroy(&rand));
 70:   }

 72:   PetscCall(VecScatterBegin(sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD));
 73:   PetscCall(VecScatterEnd(sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD));
 74:   PetscCall(VecScatterBegin(sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD));
 75:   PetscCall(VecScatterEnd(sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD));
 76:   PetscCall(VecEqual(vecLeave0, vecLeave1, &flg));
 77:   PetscCheck(flg, comm, PETSC_ERR_PLIB, "leave vectors differ");

 79:   PetscCall(VecScatterBegin(sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE));
 80:   PetscCall(VecScatterEnd(sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE));
 81:   PetscCall(VecScatterBegin(sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE));
 82:   PetscCall(VecScatterEnd(sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE));
 83:   PetscCall(VecEqual(vecRoot0, vecRoot1, &flg));
 84:   PetscCheck(flg, comm, PETSC_ERR_PLIB, "root vectors differ");

 86:   PetscCall(VecDestroy(&vecRoot0));
 87:   PetscCall(VecDestroy(&vecRoot1));
 88:   PetscCall(VecDestroy(&vecLeave0));
 89:   PetscCall(VecDestroy(&vecLeave1));
 90:   PetscFunctionReturn(PETSC_SUCCESS);
 91: }

 93: static PetscErrorCode PetscSFViewCustom(PetscSF sf, PetscViewer viewer)
 94: {
 95:   PetscMPIInt        rank, nranks;
 96:   PetscInt           i, nroots, nleaves;
 97:   const PetscInt    *ilocal;
 98:   const PetscSFNode *iremote;
 99:   PetscLayout        rootLayout;
100:   PetscInt          *gremote;

102:   PetscFunctionBegin;
103:   PetscCall(PetscSFSetUp(sf));
104:   PetscCall(PetscSFGetGraph(sf, &nroots, &nleaves, &ilocal, &iremote));
105:   PetscCall(PetscSFGetRootRanks(sf, &nranks, NULL, NULL, NULL, NULL));
106:   PetscCall(PetscSFGetGraphLayout(sf, &rootLayout, NULL, NULL, &gremote));
107:   PetscCheck(nroots == rootLayout->n, PetscObjectComm((PetscObject)sf), PETSC_ERR_PLIB, "Assertion failed: nroots == rootLayout->n");
108:   PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)sf, viewer));
109:   PetscCall(PetscViewerASCIIPushTab(viewer));
110:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)sf), &rank));
111:   PetscCall(PetscViewerASCIIPushSynchronized(viewer));
112:   if (rank == 0) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "rank #leaves #roots\n"));
113:   PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%2d] %7" PetscInt_FMT " %6" PetscInt_FMT "\n", rank, nleaves, nroots));
114:   PetscCall(PetscViewerFlush(viewer));
115:   if (rank == 0) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "leaves      roots       roots in global numbering\n"));
116:   for (i = 0; i < nleaves; i++)
117:     PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "(%2d, %2d) <- (%2d, %2" PetscInt_FMT ")  = %3" PetscInt_FMT "\n", rank, (PetscMPIInt)(ilocal ? ilocal[i] : i), (PetscMPIInt)iremote[i].rank, iremote[i].index, gremote[i]));
118:   PetscCall(PetscViewerFlush(viewer));
119:   PetscCall(PetscViewerASCIIPopSynchronized(viewer));
120:   PetscCall(PetscViewerASCIIPopTab(viewer));
121:   PetscCall(PetscLayoutDestroy(&rootLayout));
122:   PetscCall(PetscFree(gremote));
123:   PetscFunctionReturn(PETSC_SUCCESS);
124: }

126: PetscErrorCode CreateReferenceSF_Regular(AppCtx *ctx, PetscSF *refSF)
127: {
128:   PetscInt  j;
129:   PetscInt *ilocal  = NULL;
130:   PetscInt  nLeaves = ctx->nsfs * ctx->n * ctx->size;
131:   PetscInt  nroots  = ctx->n * ctx->nsfs;
132:   PetscSF   sf;

134:   PetscFunctionBegin;
135:   ilocal = NULL;
136:   if (ctx->sparseLeaves) PetscCall(PetscCalloc1(nLeaves + 1, &ilocal));
137:   PetscCall(PetscSFCreate(ctx->comm, &sf));
138:   for (j = 0; j < nLeaves; j++) {
139:     if (ctx->sparseLeaves) ilocal[j + 1] = ilocal[j] + ctx->leaveStep;
140:   }
141:   switch (ctx->rootMode) {
142:   case PETSCSF_CONCATENATE_ROOTMODE_SHARED:
143:   case PETSCSF_CONCATENATE_ROOTMODE_LOCAL: {
144:     PetscInt     i, k;
145:     PetscMPIInt  r;
146:     PetscSFNode *iremote;

148:     PetscCall(PetscCalloc1(nLeaves, &iremote));
149:     for (i = 0, j = 0; i < ctx->nsfs; i++) {
150:       for (r = 0; r < ctx->size; r++) {
151:         for (k = 0; k < ctx->n; k++, j++) {
152:           iremote[j].rank  = r;
153:           iremote[j].index = k + i * ctx->n;
154:         }
155:       }
156:     }
157:     PetscCall(PetscSFSetGraph(sf, nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER));
158:   } break;
159:   case PETSCSF_CONCATENATE_ROOTMODE_GLOBAL: {
160:     PetscLayout map = NULL;
161:     PetscInt   *gremote;

163:     PetscCall(PetscLayoutCreateFromSizes(ctx->comm, nroots, PETSC_DECIDE, 1, &map));
164:     PetscCall(PetscMalloc1(nLeaves, &gremote));
165:     for (j = 0; j < nLeaves; j++) gremote[j] = j;
166:     PetscCall(PetscSFSetGraphLayout(sf, map, nLeaves, ilocal, PETSC_OWN_POINTER, gremote));
167:     PetscCall(PetscFree(gremote));
168:     PetscCall(PetscLayoutDestroy(&map));
169:   } break;
170:   default:
171:     SETERRQ(ctx->comm, PETSC_ERR_SUP, "unsupported rootmode %d", ctx->rootMode);
172:   }
173:   PetscCall(PetscObjectSetName((PetscObject)sf, "reference_sf"));
174:   if (ctx->viewer) PetscCall(PetscSFViewCustom(sf, ctx->viewer));
175:   *refSF = sf;
176:   PetscFunctionReturn(PETSC_SUCCESS);
177: }

179: PetscErrorCode CreateSFs_Irregular(AppCtx *ctx, PetscSF *newSFs[], PetscInt *leafOffsets[])
180: {
181:   PetscInt  i;
182:   PetscInt *lOffsets = NULL;
183:   PetscSF  *sfs;
184:   PetscInt  nLeaves = ctx->n * ctx->size + (ctx->size - 1) * ctx->size / 2;
185:   PetscInt  nroots  = ctx->n + ctx->rank + ctx->nsfs - 1 + ctx->size - 1;

187:   PetscFunctionBegin;
188:   if (ctx->sparseLeaves) PetscCall(PetscCalloc1(ctx->nsfs + 1, &lOffsets));
189:   PetscCall(PetscMalloc1(ctx->nsfs, &sfs));
190:   for (i = 0; i < ctx->nsfs; i++) {
191:     PetscSF      sf;
192:     PetscInt     j, k;
193:     PetscMPIInt  r;
194:     PetscInt    *ilocal = NULL;
195:     PetscSFNode *iremote;
196:     char         name[32];

198:     if (ctx->sparseLeaves) PetscCall(PetscCalloc1(nLeaves + 1, &ilocal));
199:     PetscCall(PetscMalloc1(nLeaves, &iremote));
200:     for (r = ctx->size - 1, j = 0; r >= 0; r--) {
201:       for (k = 0; k < ctx->n + r; k++, j++) {
202:         if (ctx->sparseLeaves) ilocal[j + 1] = ilocal[j] + ctx->leaveStep;
203:         iremote[j].rank  = r;
204:         iremote[j].index = k + i + ctx->rank;
205:       }
206:     }
207:     if (ctx->sparseLeaves) lOffsets[i + 1] = lOffsets[i] + ilocal[j];

209:     PetscCall(PetscSFCreate(ctx->comm, &sf));
210:     PetscCall(PetscSFSetGraph(sf, nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER));
211:     PetscCall(PetscSNPrintf(name, sizeof(name), "sf_%" PetscInt_FMT, i));
212:     PetscCall(PetscObjectSetName((PetscObject)sf, name));
213:     if (ctx->viewer) PetscCall(PetscSFViewCustom(sf, ctx->viewer));
214:     sfs[i] = sf;
215:   }
216:   *newSFs      = sfs;
217:   *leafOffsets = lOffsets;
218:   PetscFunctionReturn(PETSC_SUCCESS);
219: }

221: PetscErrorCode CreateSFs_Regular(AppCtx *ctx, PetscSF *newSFs[], PetscInt *leafOffsets[])
222: {
223:   PetscInt                   i;
224:   PetscInt                  *lOffsets = NULL;
225:   PetscInt                   nLeaves  = ctx->n * ctx->size;
226:   PetscSF                   *sfs;
227:   PetscSFConcatenateRootMode mode = ctx->compare ? ctx->rootMode : PETSCSF_CONCATENATE_ROOTMODE_LOCAL;

229:   PetscFunctionBegin;
230:   if (ctx->sparseLeaves) PetscCall(PetscCalloc1(ctx->nsfs + 1, &lOffsets));
231:   PetscCall(PetscCalloc1(ctx->nsfs, &sfs));
232:   for (i = 0; i < ctx->nsfs; i++) {
233:     PetscSF   sf;
234:     PetscInt  j;
235:     PetscInt *ilocal = NULL;
236:     char      name[32];

238:     PetscCall(PetscSFCreate(ctx->comm, &sf));
239:     if (ctx->sparseLeaves) {
240:       PetscCall(PetscCalloc1(nLeaves + 1, &ilocal));
241:       for (j = 0; j < nLeaves; j++) ilocal[j + 1] = ilocal[j] + ctx->leaveStep;
242:       lOffsets[i + 1] = lOffsets[i] + ilocal[nLeaves];
243:     }
244:     switch (mode) {
245:     case PETSCSF_CONCATENATE_ROOTMODE_LOCAL: {
246:       PetscInt     k, nroots = ctx->n;
247:       PetscMPIInt  r;
248:       PetscSFNode *iremote;

250:       PetscCall(PetscMalloc1(nLeaves, &iremote));
251:       for (r = 0, j = 0; r < ctx->size; r++) {
252:         for (k = 0; k < ctx->n; k++, j++) {
253:           iremote[j].rank  = r;
254:           iremote[j].index = k;
255:         }
256:       }
257:       PetscCall(PetscSFSetGraph(sf, nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER));
258:     } break;
259:     case PETSCSF_CONCATENATE_ROOTMODE_SHARED: {
260:       PetscInt     k, nroots = ctx->n * ctx->nsfs;
261:       PetscMPIInt  r;
262:       PetscSFNode *iremote;

264:       PetscCall(PetscMalloc1(nLeaves, &iremote));
265:       for (r = 0, j = 0; r < ctx->size; r++) {
266:         for (k = 0; k < ctx->n; k++, j++) {
267:           iremote[j].rank  = r;
268:           iremote[j].index = k + i * ctx->n;
269:         }
270:       }
271:       PetscCall(PetscSFSetGraph(sf, nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER));
272:     } break;
273:     case PETSCSF_CONCATENATE_ROOTMODE_GLOBAL: {
274:       PetscInt    nroots = ctx->n;
275:       PetscLayout map    = NULL;
276:       PetscInt   *gremote;

278:       PetscCall(PetscLayoutCreateFromSizes(ctx->comm, nroots, PETSC_DECIDE, 1, &map));
279:       PetscCall(PetscMalloc1(nLeaves, &gremote));
280:       for (j = 0; j < nLeaves; j++) gremote[j] = j;
281:       PetscCall(PetscSFSetGraphLayout(sf, map, nLeaves, ilocal, PETSC_OWN_POINTER, gremote));
282:       PetscCall(PetscFree(gremote));
283:       PetscCall(PetscLayoutDestroy(&map));
284:     } break;
285:     default:
286:       SETERRQ(ctx->comm, PETSC_ERR_SUP, "unsupported rootmode %d", ctx->rootMode);
287:     }
288:     PetscCall(PetscSNPrintf(name, sizeof(name), "sf_%" PetscInt_FMT, i));
289:     PetscCall(PetscObjectSetName((PetscObject)sf, name));
290:     if (ctx->viewer) PetscCall(PetscSFViewCustom(sf, ctx->viewer));
291:     sfs[i] = sf;
292:   }
293:   *newSFs      = sfs;
294:   *leafOffsets = lOffsets;
295:   PetscFunctionReturn(PETSC_SUCCESS);
296: }

298: PetscErrorCode DestroySFs(AppCtx *ctx, PetscSF *sfs[])
299: {
300:   PetscInt i;

302:   PetscFunctionBegin;
303:   for (i = 0; i < ctx->nsfs; i++) PetscCall(PetscSFDestroy(&(*sfs)[i]));
304:   PetscCall(PetscFree(*sfs));
305:   PetscFunctionReturn(PETSC_SUCCESS);
306: }

308: int main(int argc, char **argv)
309: {
310:   AppCtx    ctx_;
311:   AppCtx   *ctx = &ctx_;
312:   PetscSF   sf;
313:   PetscSF  *sfs         = NULL;
314:   PetscInt *leafOffsets = NULL;
315:   MPI_Comm  comm;

317:   PetscFunctionBeginUser;
318:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
319:   comm = PETSC_COMM_WORLD;
320:   PetscCall(GetOptions(comm, ctx));

322:   if (ctx->irregular) {
323:     PetscCall(CreateSFs_Irregular(ctx, &sfs, &leafOffsets));
324:   } else {
325:     PetscCall(CreateSFs_Regular(ctx, &sfs, &leafOffsets));
326:   }
327:   PetscCall(PetscSFConcatenate(comm, ctx->nsfs, sfs, ctx->rootMode, leafOffsets, &sf));
328:   PetscCall(PetscObjectSetName((PetscObject)sf, "result_sf"));
329:   if (ctx->viewer) {
330:     PetscCall(PetscPrintf(comm, "rootMode = %s:\n", PetscSFConcatenateRootModes[ctx->rootMode]));
331:     PetscCall(PetscSFViewCustom(sf, ctx->viewer));
332:   }
333:   if (ctx->compare) {
334:     PetscSF sfRef;

336:     PetscAssert(!ctx->irregular, comm, PETSC_ERR_SUP, "Combination  -compare_to_reference true -irregular true  not implemented");
337:     PetscCall(CreateReferenceSF_Regular(ctx, &sfRef));
338:     PetscCall(PetscSFCheckEqual_Private(sf, sfRef));
339:     PetscCall(PetscSFDestroy(&sfRef));
340:   }
341:   PetscCall(DestroySFs(ctx, &sfs));
342:   PetscCall(PetscFree(leafOffsets));
343:   PetscCall(PetscSFDestroy(&sf));
344:   if (ctx->viewer) {
345:     PetscCall(PetscViewerPopFormat(ctx->viewer));
346:     PetscCall(PetscViewerDestroy(&ctx->viewer));
347:   }
348:   PetscCall(PetscFinalize());
349:   return 0;
350: }

352: /*TEST
353:   test:
354:     nsize: {{1 3}}
355:     args: -compare_to_reference -nsfs {{1 3}} -n {{0 1 5}} -leave_step {{1 3}} -root_mode {{local shared global}}

357:   test:
358:     suffix: 2
359:     nsize: 2
360:     args: -irregular {{false true}separate output} -sf_view -nsfs 3 -n 1 -leave_step {{1 3}separate output} -root_mode {{local shared global}separate output}
361: TEST*/