xref: /petsc/src/vec/is/sf/tests/ex18.c (revision 5f80ce2ab25dff0f4601e710601cbbcecf323266)
1f210f596SVaclav Hapla 
2f210f596SVaclav Hapla static char help[]= "Test PetscSFConcatenate()\n\n";
3f210f596SVaclav Hapla 
4f210f596SVaclav Hapla #include <petscsf.h>
5f210f596SVaclav Hapla 
6f210f596SVaclav Hapla typedef struct {
7f210f596SVaclav Hapla   MPI_Comm          comm;
8f210f596SVaclav Hapla   PetscMPIInt       rank, size;
9f210f596SVaclav Hapla   PetscInt          leaveStep, nsfs, nLeavesPerRank;
10f210f596SVaclav Hapla   PetscBool         shareRoots, sparseLeaves;
11f210f596SVaclav Hapla } AppCtx;
12f210f596SVaclav Hapla 
13f210f596SVaclav Hapla static PetscErrorCode GetOptions(MPI_Comm comm, AppCtx *ctx)
14f210f596SVaclav Hapla {
15f210f596SVaclav Hapla   PetscFunctionBegin;
16f210f596SVaclav Hapla   ctx->comm = comm;
17f210f596SVaclav Hapla   ctx->nsfs = 3;
18f210f596SVaclav Hapla   ctx->nLeavesPerRank = 4;
19f210f596SVaclav Hapla   ctx->leaveStep = 1;
20f210f596SVaclav Hapla   ctx->shareRoots = PETSC_FALSE;
21f210f596SVaclav Hapla   ctx->sparseLeaves = PETSC_FALSE;
22*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsGetInt(NULL, NULL, "-nsfs", &ctx->nsfs, NULL));
23*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsGetInt(NULL, NULL, "-n_leaves_per_rank", &ctx->nLeavesPerRank, NULL));
24*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsGetInt(NULL, NULL, "-leave_step", &ctx->leaveStep, NULL));
25*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsGetBool(NULL, NULL, "-share_roots", &ctx->shareRoots, NULL));
26f210f596SVaclav Hapla   ctx->sparseLeaves = (PetscBool) (ctx->leaveStep != 1);
27*5f80ce2aSJacob Faibussowitsch   CHKERRMPI(MPI_Comm_size(comm, &ctx->size));
28*5f80ce2aSJacob Faibussowitsch   CHKERRMPI(MPI_Comm_rank(comm, &ctx->rank));
29f210f596SVaclav Hapla   PetscFunctionReturn(0);
30f210f596SVaclav Hapla }
31f210f596SVaclav Hapla 
32f210f596SVaclav Hapla static PetscErrorCode PetscSFCheckEqual_Private(PetscSF sf0, PetscSF sf1)
33f210f596SVaclav Hapla {
34f210f596SVaclav Hapla   PetscInt          nRoot, nLeave;
35f210f596SVaclav Hapla   Vec               vecRoot0, vecLeave0, vecRoot1, vecLeave1;
36f210f596SVaclav Hapla   MPI_Comm          comm;
37f210f596SVaclav Hapla   PetscBool         flg;
38f210f596SVaclav Hapla 
39f210f596SVaclav Hapla   PetscFunctionBegin;
40*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscObjectGetComm((PetscObject)sf0, &comm));
41*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscSFGetGraph(sf0, &nRoot, NULL, NULL, NULL));
42*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscSFGetLeafRange(sf0, NULL, &nLeave));
43f210f596SVaclav Hapla   nLeave++;
44*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecCreateMPI(comm, nRoot, PETSC_DECIDE, &vecRoot0));
45*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecCreateMPI(comm, nLeave, PETSC_DECIDE, &vecLeave0));
46*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecDuplicate(vecRoot0, &vecRoot1));
47*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecDuplicate(vecLeave0, &vecLeave1));
48f210f596SVaclav Hapla   {
49f210f596SVaclav Hapla     PetscRandom       rand;
50f210f596SVaclav Hapla 
51*5f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscRandomCreate(comm, &rand));
52*5f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscRandomSetFromOptions(rand));
53*5f80ce2aSJacob Faibussowitsch     CHKERRQ(VecSetRandom(vecRoot0, rand));
54*5f80ce2aSJacob Faibussowitsch     CHKERRQ(VecSetRandom(vecLeave0, rand));
55*5f80ce2aSJacob Faibussowitsch     CHKERRQ(VecCopy(vecRoot0, vecRoot1));
56*5f80ce2aSJacob Faibussowitsch     CHKERRQ(VecCopy(vecLeave0, vecLeave1));
57*5f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscRandomDestroy(&rand));
58f210f596SVaclav Hapla   }
59f210f596SVaclav Hapla 
60*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecScatterBegin(sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD));
61*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecScatterEnd(  sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD));
62*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecScatterBegin(sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD));
63*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecScatterEnd(  sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD));
64*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecEqual(vecLeave0, vecLeave1, &flg));
65f210f596SVaclav Hapla   PetscCheck(flg, comm, PETSC_ERR_PLIB, "leave vectors differ");
66f210f596SVaclav Hapla 
67*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecScatterBegin(sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE));
68*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecScatterEnd(  sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE));
69*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecScatterBegin(sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE));
70*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecScatterEnd(  sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE));
71*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecEqual(vecRoot0, vecRoot1, &flg));
72f210f596SVaclav Hapla   PetscCheck(flg, comm, PETSC_ERR_PLIB, "root vectors differ");
73f210f596SVaclav Hapla 
74*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecDestroy(&vecRoot0));
75*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecDestroy(&vecRoot1));
76*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecDestroy(&vecLeave0));
77*5f80ce2aSJacob Faibussowitsch   CHKERRQ(VecDestroy(&vecLeave1));
78f210f596SVaclav Hapla   PetscFunctionReturn(0);
79f210f596SVaclav Hapla }
80f210f596SVaclav Hapla 
81f210f596SVaclav Hapla PetscErrorCode CreateReferenceSF(AppCtx *ctx, PetscSF *refSF)
82f210f596SVaclav Hapla {
83f210f596SVaclav Hapla   PetscInt          i, j, k, r;
84f210f596SVaclav Hapla   PetscInt         *ilocal = NULL;
85f210f596SVaclav Hapla   PetscSFNode      *iremote;
86f210f596SVaclav Hapla   PetscInt          nLeaves = ctx->nsfs * ctx->nLeavesPerRank * ctx->size;
87f210f596SVaclav Hapla   PetscInt          nroots  = ctx->nLeavesPerRank * ctx->nsfs;
88f210f596SVaclav Hapla   PetscSF           sf;
89f210f596SVaclav Hapla 
90f210f596SVaclav Hapla   PetscFunctionBegin;
91f210f596SVaclav Hapla   ilocal = NULL;
92f210f596SVaclav Hapla   if (ctx->sparseLeaves) {
93*5f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscCalloc1(nLeaves+1, &ilocal));
94f210f596SVaclav Hapla   }
95*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscMalloc1(nLeaves, &iremote));
96*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscSFCreate(ctx->comm, &sf));
97f210f596SVaclav Hapla   for (i=0, j=0; i<ctx->nsfs; i++) {
98f210f596SVaclav Hapla     for (r=0; r<ctx->size; r++) {
99f210f596SVaclav Hapla       for (k=0; k<ctx->nLeavesPerRank; k++, j++) {
100f210f596SVaclav Hapla         if (ctx->sparseLeaves) {
101f210f596SVaclav Hapla           ilocal[j+1] = ilocal[j] + ctx->leaveStep;
102f210f596SVaclav Hapla         }
103f210f596SVaclav Hapla         iremote[j].rank = r;
104f210f596SVaclav Hapla         iremote[j].index = k + i * ctx->nLeavesPerRank;
105f210f596SVaclav Hapla       }
106f210f596SVaclav Hapla     }
107f210f596SVaclav Hapla   }
108*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscSFSetGraph(sf, nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER));
109f210f596SVaclav Hapla   *refSF = sf;
110f210f596SVaclav Hapla   PetscFunctionReturn(0);
111f210f596SVaclav Hapla }
112f210f596SVaclav Hapla 
113f210f596SVaclav Hapla PetscErrorCode CreateSFs(AppCtx *ctx, PetscSF *newSFs[], PetscInt *leafOffsets[])
114f210f596SVaclav Hapla {
115f210f596SVaclav Hapla   PetscInt          i;
116f210f596SVaclav Hapla   PetscInt         *lOffsets = NULL;
117f210f596SVaclav Hapla   PetscSF          *sfs;
118f210f596SVaclav Hapla   PetscInt          nLeaves = ctx->nLeavesPerRank * ctx->size;
119f210f596SVaclav Hapla   PetscInt          nroots  = ctx->shareRoots ? ctx->nLeavesPerRank * ctx->nsfs : ctx->nLeavesPerRank;
120f210f596SVaclav Hapla 
121f210f596SVaclav Hapla   PetscFunctionBegin;
122f210f596SVaclav Hapla   if (ctx->sparseLeaves) {
123*5f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscCalloc1(ctx->nsfs+1, &lOffsets));
124f210f596SVaclav Hapla   }
125*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscMalloc1(ctx->nsfs, &sfs));
126f210f596SVaclav Hapla   for (i=0; i<ctx->nsfs; i++) {
127f210f596SVaclav Hapla     PetscInt      j, k;
128f210f596SVaclav Hapla     PetscMPIInt   r;
129f210f596SVaclav Hapla     PetscInt     *ilocal = NULL;
130f210f596SVaclav Hapla     PetscSFNode  *iremote;
131f210f596SVaclav Hapla 
132f210f596SVaclav Hapla     if (ctx->sparseLeaves) {
133*5f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscCalloc1(nLeaves+1, &ilocal));
134f210f596SVaclav Hapla     }
135*5f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscMalloc1(nLeaves, &iremote));
136f210f596SVaclav Hapla     for (r=0, j=0; r<ctx->size; r++) {
137f210f596SVaclav Hapla       for (k=0; k<ctx->nLeavesPerRank; k++, j++) {
138f210f596SVaclav Hapla         if (ctx->sparseLeaves) {
139f210f596SVaclav Hapla           ilocal[j+1] = ilocal[j] + ctx->leaveStep;
140f210f596SVaclav Hapla         }
141f210f596SVaclav Hapla         iremote[j].rank = r;
142f210f596SVaclav Hapla         iremote[j].index = ctx->shareRoots ? k + i * ctx->nLeavesPerRank : k;
143f210f596SVaclav Hapla       }
144f210f596SVaclav Hapla     }
145f210f596SVaclav Hapla     if (ctx->sparseLeaves) lOffsets[i+1] = lOffsets[i] + ilocal[j];
146f210f596SVaclav Hapla 
147*5f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscSFCreate(ctx->comm, &sfs[i]));
148*5f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscSFSetGraph(sfs[i], nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER));
149f210f596SVaclav Hapla   }
150f210f596SVaclav Hapla   *newSFs = sfs;
151f210f596SVaclav Hapla   *leafOffsets = lOffsets;
152f210f596SVaclav Hapla   PetscFunctionReturn(0);
153f210f596SVaclav Hapla }
154f210f596SVaclav Hapla 
155f210f596SVaclav Hapla PetscErrorCode DestroySFs(AppCtx *ctx, PetscSF *sfs[])
156f210f596SVaclav Hapla {
157f210f596SVaclav Hapla   PetscInt          i;
158f210f596SVaclav Hapla 
159f210f596SVaclav Hapla   PetscFunctionBegin;
160f210f596SVaclav Hapla   for (i=0; i<ctx->nsfs; i++) {
161*5f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscSFDestroy(&(*sfs)[i]));
162f210f596SVaclav Hapla   }
163*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscFree(*sfs));
164f210f596SVaclav Hapla   PetscFunctionReturn(0);
165f210f596SVaclav Hapla }
166f210f596SVaclav Hapla 
167f210f596SVaclav Hapla int main(int argc, char **argv)
168f210f596SVaclav Hapla {
169f210f596SVaclav Hapla   AppCtx            ctx;
170f210f596SVaclav Hapla   PetscSF           sf, sfRef;
171f210f596SVaclav Hapla   PetscSF          *sfs;
172f210f596SVaclav Hapla   PetscInt         *leafOffsets;
173f210f596SVaclav Hapla   MPI_Comm          comm;
174f210f596SVaclav Hapla   PetscErrorCode    ierr;
175f210f596SVaclav Hapla 
176f210f596SVaclav Hapla   ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr;
177f210f596SVaclav Hapla   comm = PETSC_COMM_WORLD;
178*5f80ce2aSJacob Faibussowitsch   CHKERRQ(GetOptions(comm, &ctx));
179f210f596SVaclav Hapla 
180*5f80ce2aSJacob Faibussowitsch   CHKERRQ(CreateSFs(&ctx, &sfs, &leafOffsets));
181*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscSFConcatenate(comm, ctx.nsfs, sfs, ctx.shareRoots, leafOffsets, &sf));
182*5f80ce2aSJacob Faibussowitsch   CHKERRQ(CreateReferenceSF(&ctx, &sfRef));
183*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscSFCheckEqual_Private(sf, sfRef));
184f210f596SVaclav Hapla 
185*5f80ce2aSJacob Faibussowitsch   CHKERRQ(DestroySFs(&ctx, &sfs));
186*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscFree(leafOffsets));
187*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscSFDestroy(&sf));
188*5f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscSFDestroy(&sfRef));
189f210f596SVaclav Hapla   ierr = PetscFinalize();
190f210f596SVaclav Hapla   return ierr;
191f210f596SVaclav Hapla }
192f210f596SVaclav Hapla 
193f210f596SVaclav Hapla /*TEST
194f210f596SVaclav Hapla   test:
195f210f596SVaclav Hapla     nsize: {{1 3}}
196f210f596SVaclav Hapla     args: -nsfs {{1 3}} -n_leaves_per_rank {{0 1 5}} -leave_step {{1 3}} -share_roots {{true false}}
197f210f596SVaclav Hapla TEST*/
198