1 static char help[] = "Test PetscSFConcatenate()\n\n";
2
3 #include <petscsf.h>
4
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;
15
GetOptions(MPI_Comm comm,AppCtx * ctx)16 static PetscErrorCode GetOptions(MPI_Comm comm, AppCtx *ctx)
17 {
18 PetscViewerFormat format;
19
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 }
43
PetscSFCheckEqual_Private(PetscSF sf0,PetscSF sf1)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;
50
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;
62
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 }
71
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");
78
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");
85
86 PetscCall(VecDestroy(&vecRoot0));
87 PetscCall(VecDestroy(&vecRoot1));
88 PetscCall(VecDestroy(&vecLeave0));
89 PetscCall(VecDestroy(&vecLeave1));
90 PetscFunctionReturn(PETSC_SUCCESS);
91 }
92
PetscSFViewCustom(PetscSF sf,PetscViewer viewer)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;
101
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 }
125
CreateReferenceSF_Regular(AppCtx * ctx,PetscSF * refSF)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;
133
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;
147
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;
162
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 }
178
CreateSFs_Irregular(AppCtx * ctx,PetscSF * newSFs[],PetscInt * leafOffsets[])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;
186
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];
197
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];
208
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 }
220
CreateSFs_Regular(AppCtx * ctx,PetscSF * newSFs[],PetscInt * leafOffsets[])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;
228
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];
237
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;
249
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;
263
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;
277
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 }
297
DestroySFs(AppCtx * ctx,PetscSF * sfs[])298 PetscErrorCode DestroySFs(AppCtx *ctx, PetscSF *sfs[])
299 {
300 PetscInt i;
301
302 PetscFunctionBegin;
303 for (i = 0; i < ctx->nsfs; i++) PetscCall(PetscSFDestroy(&(*sfs)[i]));
304 PetscCall(PetscFree(*sfs));
305 PetscFunctionReturn(PETSC_SUCCESS);
306 }
307
main(int argc,char ** argv)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;
316
317 PetscFunctionBeginUser;
318 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
319 comm = PETSC_COMM_WORLD;
320 PetscCall(GetOptions(comm, ctx));
321
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;
335
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 }
351
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}}
356 output_file: output/empty.out
357
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*/
363