Actual source code: ex18.c


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

  4: #include <petscsf.h>

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

358:   test:
359:     suffix: 2
360:     nsize: 2
361:     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}
362: TEST*/