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) {
90: PetscCalloc1(nLeaves+1, &ilocal);
91: }
92: PetscMalloc1(nLeaves, &iremote);
93: PetscSFCreate(ctx->comm, &sf);
94: for (i=0, j=0; i<ctx->nsfs; i++) {
95: for (r=0; r<ctx->size; r++) {
96: for (k=0; k<ctx->nLeavesPerRank; k++, j++) {
97: if (ctx->sparseLeaves) {
98: ilocal[j+1] = ilocal[j] + ctx->leaveStep;
99: }
100: iremote[j].rank = r;
101: iremote[j].index = k + i * ctx->nLeavesPerRank;
102: }
103: }
104: }
105: PetscSFSetGraph(sf, nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER);
106: *refSF = sf;
107: return 0;
108: }
110: PetscErrorCode CreateSFs(AppCtx *ctx, PetscSF *newSFs[], PetscInt *leafOffsets[])
111: {
112: PetscInt i;
113: PetscInt *lOffsets = NULL;
114: PetscSF *sfs;
115: PetscInt nLeaves = ctx->nLeavesPerRank * ctx->size;
116: PetscInt nroots = ctx->shareRoots ? ctx->nLeavesPerRank * ctx->nsfs : ctx->nLeavesPerRank;
118: if (ctx->sparseLeaves) {
119: PetscCalloc1(ctx->nsfs+1, &lOffsets);
120: }
121: PetscMalloc1(ctx->nsfs, &sfs);
122: for (i=0; i<ctx->nsfs; i++) {
123: PetscInt j, k;
124: PetscMPIInt r;
125: PetscInt *ilocal = NULL;
126: PetscSFNode *iremote;
128: if (ctx->sparseLeaves) {
129: PetscCalloc1(nLeaves+1, &ilocal);
130: }
131: PetscMalloc1(nLeaves, &iremote);
132: for (r=0, j=0; r<ctx->size; r++) {
133: for (k=0; k<ctx->nLeavesPerRank; k++, j++) {
134: if (ctx->sparseLeaves) {
135: ilocal[j+1] = ilocal[j] + ctx->leaveStep;
136: }
137: iremote[j].rank = r;
138: iremote[j].index = ctx->shareRoots ? k + i * ctx->nLeavesPerRank : k;
139: }
140: }
141: if (ctx->sparseLeaves) lOffsets[i+1] = lOffsets[i] + ilocal[j];
143: PetscSFCreate(ctx->comm, &sfs[i]);
144: PetscSFSetGraph(sfs[i], nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER);
145: }
146: *newSFs = sfs;
147: *leafOffsets = lOffsets;
148: return 0;
149: }
151: PetscErrorCode DestroySFs(AppCtx *ctx, PetscSF *sfs[])
152: {
153: PetscInt i;
155: for (i=0; i<ctx->nsfs; i++) {
156: PetscSFDestroy(&(*sfs)[i]);
157: }
158: PetscFree(*sfs);
159: return 0;
160: }
162: int main(int argc, char **argv)
163: {
164: AppCtx ctx;
165: PetscSF sf, sfRef;
166: PetscSF *sfs;
167: PetscInt *leafOffsets;
168: MPI_Comm comm;
170: PetscInitialize(&argc,&argv,NULL,help);
171: comm = PETSC_COMM_WORLD;
172: GetOptions(comm, &ctx);
174: CreateSFs(&ctx, &sfs, &leafOffsets);
175: PetscSFConcatenate(comm, ctx.nsfs, sfs, ctx.shareRoots, leafOffsets, &sf);
176: CreateReferenceSF(&ctx, &sfRef);
177: PetscSFCheckEqual_Private(sf, sfRef);
179: DestroySFs(&ctx, &sfs);
180: PetscFree(leafOffsets);
181: PetscSFDestroy(&sf);
182: PetscSFDestroy(&sfRef);
183: PetscFinalize();
184: return 0;
185: }
187: /*TEST
188: test:
189: nsize: {{1 3}}
190: args: -nsfs {{1 3}} -n_leaves_per_rank {{0 1 5}} -leave_step {{1 3}} -share_roots {{true false}}
191: TEST*/