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