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, nLeavesPerRank;
10: PetscBool shareRoots, sparseLeaves;
11: } AppCtx;
13: static PetscErrorCode GetOptions(MPI_Comm comm, AppCtx *ctx)
14: {
15: ctx->comm = comm;
16: ctx->nsfs = 3;
17: ctx->nLeavesPerRank = 4;
18: ctx->leaveStep = 1;
19: ctx->shareRoots = PETSC_FALSE;
20: ctx->sparseLeaves = PETSC_FALSE;
21: PetscOptionsGetInt(NULL, NULL, "-nsfs", &ctx->nsfs, NULL);
22: PetscOptionsGetInt(NULL, NULL, "-n_leaves_per_rank", &ctx->nLeavesPerRank, NULL);
23: PetscOptionsGetInt(NULL, NULL, "-leave_step", &ctx->leaveStep, NULL);
24: PetscOptionsGetBool(NULL, NULL, "-share_roots", &ctx->shareRoots, NULL);
25: ctx->sparseLeaves = (PetscBool)(ctx->leaveStep != 1);
26: MPI_Comm_size(comm, &ctx->size);
27: MPI_Comm_rank(comm, &ctx->rank);
28: return 0;
29: }
31: static PetscErrorCode PetscSFCheckEqual_Private(PetscSF sf0, PetscSF sf1)
32: {
33: PetscInt nRoot, nLeave;
34: Vec vecRoot0, vecLeave0, vecRoot1, vecLeave1;
35: MPI_Comm comm;
36: PetscBool flg;
38: PetscObjectGetComm((PetscObject)sf0, &comm);
39: PetscSFGetGraph(sf0, &nRoot, NULL, NULL, NULL);
40: PetscSFGetLeafRange(sf0, NULL, &nLeave);
41: nLeave++;
42: VecCreateMPI(comm, nRoot, PETSC_DECIDE, &vecRoot0);
43: VecCreateMPI(comm, nLeave, PETSC_DECIDE, &vecLeave0);
44: VecDuplicate(vecRoot0, &vecRoot1);
45: VecDuplicate(vecLeave0, &vecLeave1);
46: {
47: PetscRandom rand;
49: PetscRandomCreate(comm, &rand);
50: PetscRandomSetFromOptions(rand);
51: VecSetRandom(vecRoot0, rand);
52: VecSetRandom(vecLeave0, rand);
53: VecCopy(vecRoot0, vecRoot1);
54: VecCopy(vecLeave0, vecLeave1);
55: PetscRandomDestroy(&rand);
56: }
58: VecScatterBegin(sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD);
59: VecScatterEnd(sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD);
60: VecScatterBegin(sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD);
61: VecScatterEnd(sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD);
62: VecEqual(vecLeave0, vecLeave1, &flg);
65: VecScatterBegin(sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE);
66: VecScatterEnd(sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE);
67: VecScatterBegin(sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE);
68: VecScatterEnd(sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE);
69: VecEqual(vecRoot0, vecRoot1, &flg);
72: VecDestroy(&vecRoot0);
73: VecDestroy(&vecRoot1);
74: VecDestroy(&vecLeave0);
75: VecDestroy(&vecLeave1);
76: return 0;
77: }
79: PetscErrorCode CreateReferenceSF(AppCtx *ctx, PetscSF *refSF)
80: {
81: PetscInt i, j, k, r;
82: PetscInt *ilocal = NULL;
83: PetscSFNode *iremote;
84: PetscInt nLeaves = ctx->nsfs * ctx->nLeavesPerRank * ctx->size;
85: PetscInt nroots = ctx->nLeavesPerRank * ctx->nsfs;
86: PetscSF sf;
88: ilocal = NULL;
89: if (ctx->sparseLeaves) PetscCalloc1(nLeaves + 1, &ilocal);
90: PetscMalloc1(nLeaves, &iremote);
91: PetscSFCreate(ctx->comm, &sf);
92: for (i = 0, j = 0; i < ctx->nsfs; i++) {
93: for (r = 0; r < ctx->size; r++) {
94: for (k = 0; k < ctx->nLeavesPerRank; k++, j++) {
95: if (ctx->sparseLeaves) ilocal[j + 1] = ilocal[j] + ctx->leaveStep;
96: iremote[j].rank = r;
97: iremote[j].index = k + i * ctx->nLeavesPerRank;
98: }
99: }
100: }
101: PetscSFSetGraph(sf, nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER);
102: *refSF = sf;
103: return 0;
104: }
106: PetscErrorCode CreateSFs(AppCtx *ctx, PetscSF *newSFs[], PetscInt *leafOffsets[])
107: {
108: PetscInt i;
109: PetscInt *lOffsets = NULL;
110: PetscSF *sfs;
111: PetscInt nLeaves = ctx->nLeavesPerRank * ctx->size;
112: PetscInt nroots = ctx->shareRoots ? ctx->nLeavesPerRank * ctx->nsfs : ctx->nLeavesPerRank;
114: if (ctx->sparseLeaves) PetscCalloc1(ctx->nsfs + 1, &lOffsets);
115: PetscMalloc1(ctx->nsfs, &sfs);
116: for (i = 0; i < ctx->nsfs; i++) {
117: PetscInt j, k;
118: PetscMPIInt r;
119: PetscInt *ilocal = NULL;
120: PetscSFNode *iremote;
122: if (ctx->sparseLeaves) PetscCalloc1(nLeaves + 1, &ilocal);
123: PetscMalloc1(nLeaves, &iremote);
124: for (r = 0, j = 0; r < ctx->size; r++) {
125: for (k = 0; k < ctx->nLeavesPerRank; k++, j++) {
126: if (ctx->sparseLeaves) ilocal[j + 1] = ilocal[j] + ctx->leaveStep;
127: iremote[j].rank = r;
128: iremote[j].index = ctx->shareRoots ? k + i * ctx->nLeavesPerRank : k;
129: }
130: }
131: if (ctx->sparseLeaves) lOffsets[i + 1] = lOffsets[i] + ilocal[j];
133: PetscSFCreate(ctx->comm, &sfs[i]);
134: PetscSFSetGraph(sfs[i], nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER);
135: }
136: *newSFs = sfs;
137: *leafOffsets = lOffsets;
138: return 0;
139: }
141: PetscErrorCode DestroySFs(AppCtx *ctx, PetscSF *sfs[])
142: {
143: PetscInt i;
145: for (i = 0; i < ctx->nsfs; i++) PetscSFDestroy(&(*sfs)[i]);
146: PetscFree(*sfs);
147: return 0;
148: }
150: int main(int argc, char **argv)
151: {
152: AppCtx ctx;
153: PetscSF sf, sfRef;
154: PetscSF *sfs;
155: PetscInt *leafOffsets;
156: MPI_Comm comm;
159: PetscInitialize(&argc, &argv, NULL, help);
160: comm = PETSC_COMM_WORLD;
161: GetOptions(comm, &ctx);
163: CreateSFs(&ctx, &sfs, &leafOffsets);
164: PetscSFConcatenate(comm, ctx.nsfs, sfs, ctx.shareRoots, leafOffsets, &sf);
165: CreateReferenceSF(&ctx, &sfRef);
166: PetscSFCheckEqual_Private(sf, sfRef);
168: DestroySFs(&ctx, &sfs);
169: PetscFree(leafOffsets);
170: PetscSFDestroy(&sf);
171: PetscSFDestroy(&sfRef);
172: PetscFinalize();
173: return 0;
174: }
176: /*TEST
177: test:
178: nsize: {{1 3}}
179: args: -nsfs {{1 3}} -n_leaves_per_rank {{0 1 5}} -leave_step {{1 3}} -share_roots {{true false}}
180: TEST*/