1 static char help[] = "Tests for DMLabel\n\n";
2
3 #include <petscdmplex.h>
4 #include <petsc/private/dmimpl.h>
5
TestInsertion()6 static PetscErrorCode TestInsertion()
7 {
8 DMLabel label, label2;
9 const PetscInt values[5] = {0, 3, 4, -1, 176}, N = 10000;
10 PetscInt i, v;
11
12 PetscFunctionBegin;
13 PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Test Label", &label));
14 PetscCall(DMLabelSetDefaultValue(label, -100));
15 for (i = 0; i < N; ++i) PetscCall(DMLabelSetValue(label, i, values[i % 5]));
16 /* Test get in hash mode */
17 for (i = 0; i < N; ++i) {
18 PetscInt val;
19
20 PetscCall(DMLabelGetValue(label, i, &val));
21 PetscCheck(val == values[i % 5], PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Value %" PetscInt_FMT " for point %" PetscInt_FMT " should be %" PetscInt_FMT, val, i, values[i % 5]);
22 }
23 /* Test stratum */
24 for (v = 0; v < 5; ++v) {
25 IS stratum;
26 const PetscInt *points;
27 PetscInt n;
28
29 PetscCall(DMLabelGetStratumIS(label, values[v], &stratum));
30 PetscCheck(stratum, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Stratum %" PetscInt_FMT " is empty!", v);
31 PetscCall(ISGetIndices(stratum, &points));
32 PetscCall(ISGetLocalSize(stratum, &n));
33 for (i = 0; i < n; ++i) PetscCheck(points[i] == i * 5 + v, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Point %" PetscInt_FMT " should be %" PetscInt_FMT, points[i], i * 5 + v);
34 PetscCall(ISRestoreIndices(stratum, &points));
35 PetscCall(ISDestroy(&stratum));
36 }
37 /* Test get in array mode */
38 for (i = 0; i < N; ++i) {
39 PetscInt val;
40
41 PetscCall(DMLabelGetValue(label, i, &val));
42 PetscCheck(val == values[i % 5], PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Value %" PetscInt_FMT " should be %" PetscInt_FMT, val, values[i % 5]);
43 }
44 /* Test Duplicate */
45 PetscCall(DMLabelDuplicate(label, &label2));
46 for (i = 0; i < N; ++i) {
47 PetscInt val;
48
49 PetscCall(DMLabelGetValue(label2, i, &val));
50 PetscCheck(val == values[i % 5], PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Value %" PetscInt_FMT " should be %" PetscInt_FMT, val, values[i % 5]);
51 }
52 PetscCall(DMLabelDestroy(&label2));
53 PetscCall(DMLabelDestroy(&label));
54 PetscFunctionReturn(PETSC_SUCCESS);
55 }
56
TestEmptyStrata(MPI_Comm comm)57 static PetscErrorCode TestEmptyStrata(MPI_Comm comm)
58 {
59 DM dm, dmDist;
60 PetscPartitioner part;
61 PetscInt c0[6] = {2, 3, 6, 7, 9, 11};
62 PetscInt c1[6] = {4, 5, 7, 8, 10, 12};
63 PetscInt c2[4] = {13, 15, 19, 21};
64 PetscInt c3[4] = {14, 16, 20, 22};
65 PetscInt c4[4] = {15, 17, 21, 23};
66 PetscInt c5[4] = {16, 18, 22, 24};
67 PetscInt c6[4] = {13, 14, 19, 20};
68 PetscInt c7[4] = {15, 16, 21, 22};
69 PetscInt c8[4] = {17, 18, 23, 24};
70 PetscInt c9[4] = {13, 14, 15, 16};
71 PetscInt c10[4] = {15, 16, 17, 18};
72 PetscInt c11[4] = {19, 20, 21, 22};
73 PetscInt c12[4] = {21, 22, 23, 24};
74 PetscInt dim = 3;
75 PetscMPIInt rank;
76
77 PetscFunctionBegin;
78 PetscCallMPI(MPI_Comm_rank(comm, &rank));
79 /* A 3D box with two adjacent cells, sharing one face and four vertices */
80 PetscCall(DMCreate(comm, &dm));
81 PetscCall(DMSetType(dm, DMPLEX));
82 PetscCall(DMSetDimension(dm, dim));
83 if (rank == 0) {
84 PetscCall(DMPlexSetChart(dm, 0, 25));
85 PetscCall(DMPlexSetConeSize(dm, 0, 6));
86 PetscCall(DMPlexSetConeSize(dm, 1, 6));
87 PetscCall(DMPlexSetConeSize(dm, 2, 4));
88 PetscCall(DMPlexSetConeSize(dm, 3, 4));
89 PetscCall(DMPlexSetConeSize(dm, 4, 4));
90 PetscCall(DMPlexSetConeSize(dm, 5, 4));
91 PetscCall(DMPlexSetConeSize(dm, 6, 4));
92 PetscCall(DMPlexSetConeSize(dm, 7, 4));
93 PetscCall(DMPlexSetConeSize(dm, 8, 4));
94 PetscCall(DMPlexSetConeSize(dm, 9, 4));
95 PetscCall(DMPlexSetConeSize(dm, 10, 4));
96 PetscCall(DMPlexSetConeSize(dm, 11, 4));
97 PetscCall(DMPlexSetConeSize(dm, 12, 4));
98 }
99 PetscCall(DMSetUp(dm));
100 if (rank == 0) {
101 PetscCall(DMPlexSetCone(dm, 0, c0));
102 PetscCall(DMPlexSetCone(dm, 1, c1));
103 PetscCall(DMPlexSetCone(dm, 2, c2));
104 PetscCall(DMPlexSetCone(dm, 3, c3));
105 PetscCall(DMPlexSetCone(dm, 4, c4));
106 PetscCall(DMPlexSetCone(dm, 5, c5));
107 PetscCall(DMPlexSetCone(dm, 6, c6));
108 PetscCall(DMPlexSetCone(dm, 7, c7));
109 PetscCall(DMPlexSetCone(dm, 8, c8));
110 PetscCall(DMPlexSetCone(dm, 9, c9));
111 PetscCall(DMPlexSetCone(dm, 10, c10));
112 PetscCall(DMPlexSetCone(dm, 11, c11));
113 PetscCall(DMPlexSetCone(dm, 12, c12));
114 }
115 PetscCall(DMPlexSymmetrize(dm));
116 /* Create a user managed depth label, so that we can leave out edges */
117 {
118 DMLabel label;
119 PetscInt numValues, maxValues = 0, v;
120
121 PetscCall(DMCreateLabel(dm, "depth"));
122 PetscCall(DMPlexGetDepthLabel(dm, &label));
123 if (rank == 0) {
124 PetscInt i;
125
126 for (i = 0; i < 25; ++i) {
127 if (i < 2) PetscCall(DMLabelSetValue(label, i, 3));
128 else if (i < 13) PetscCall(DMLabelSetValue(label, i, 2));
129 else {
130 if (i == 13) PetscCall(DMLabelAddStratum(label, 1));
131 PetscCall(DMLabelSetValue(label, i, 0));
132 }
133 }
134 }
135 PetscCall(DMLabelGetNumValues(label, &numValues));
136 PetscCallMPI(MPIU_Allreduce(&numValues, &maxValues, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)dm)));
137 for (v = numValues; v < maxValues; ++v) PetscCall(DMLabelAddStratum(label, v));
138 }
139 {
140 DMLabel label;
141 PetscCall(DMPlexGetDepthLabel(dm, &label));
142 PetscCall(DMLabelView(label, PETSC_VIEWER_STDOUT_(comm)));
143 }
144 PetscCall(DMPlexGetPartitioner(dm, &part));
145 PetscCall(PetscPartitionerSetFromOptions(part));
146 PetscCall(DMPlexDistribute(dm, 1, NULL, &dmDist));
147 if (dmDist) {
148 PetscCall(DMDestroy(&dm));
149 dm = dmDist;
150 }
151 {
152 DMLabel label;
153 PetscCall(DMPlexGetDepthLabel(dm, &label));
154 PetscCall(DMLabelView(label, PETSC_VIEWER_STDOUT_(comm)));
155 }
156 /* Create a cell vector */
157 {
158 Vec v;
159 PetscSection s;
160 PetscInt numComp[] = {1};
161 PetscInt dof[] = {0, 0, 0, 1};
162 PetscInt N;
163
164 PetscCall(DMSetNumFields(dm, 1));
165 PetscCall(DMPlexCreateSection(dm, NULL, numComp, dof, 0, NULL, NULL, NULL, NULL, &s));
166 PetscCall(DMSetLocalSection(dm, s));
167 PetscCall(PetscSectionDestroy(&s));
168 PetscCall(DMCreateGlobalVector(dm, &v));
169 PetscCall(VecGetSize(v, &N));
170 if (N != 2) {
171 PetscCall(DMView(dm, PETSC_VIEWER_STDOUT_(comm)));
172 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "FAIL: Vector size %" PetscInt_FMT " != 2", N);
173 }
174 PetscCall(VecDestroy(&v));
175 }
176 PetscCall(DMDestroy(&dm));
177 PetscFunctionReturn(PETSC_SUCCESS);
178 }
179
TestDistribution(MPI_Comm comm)180 static PetscErrorCode TestDistribution(MPI_Comm comm)
181 {
182 DM dm, dmDist;
183 PetscPartitioner part;
184 DMLabel label;
185 char filename[PETSC_MAX_PATH_LEN];
186 const char *name = "test label";
187 PetscInt overlap = 0, cStart, cEnd, c;
188 PetscMPIInt rank;
189 PetscBool flg;
190
191 PetscFunctionBegin;
192 PetscCallMPI(MPI_Comm_rank(comm, &rank));
193 PetscCall(PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), &flg));
194 if (!flg) PetscFunctionReturn(PETSC_SUCCESS);
195 PetscCall(PetscOptionsGetInt(NULL, NULL, "-overlap", &overlap, NULL));
196 PetscCall(DMPlexCreateFromFile(comm, filename, "ex11_plex", PETSC_TRUE, &dm));
197 PetscCall(DMSetBasicAdjacency(dm, PETSC_TRUE, PETSC_FALSE));
198 PetscCall(DMCreateLabel(dm, name));
199 PetscCall(DMGetLabel(dm, name, &label));
200 PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd));
201 for (c = cStart; c < cEnd; ++c) PetscCall(DMLabelSetValue(label, c, c));
202 PetscCall(DMLabelView(label, PETSC_VIEWER_STDOUT_WORLD));
203 PetscCall(DMPlexGetPartitioner(dm, &part));
204 PetscCall(PetscPartitionerSetFromOptions(part));
205 PetscCall(DMPlexDistribute(dm, overlap, NULL, &dmDist));
206 if (dmDist) {
207 PetscCall(DMDestroy(&dm));
208 dm = dmDist;
209 }
210 PetscCall(PetscObjectSetName((PetscObject)dm, "Mesh"));
211 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
212 PetscCall(DMGetLabel(dm, name, &label));
213 PetscCall(DMLabelView(label, PETSC_VIEWER_STDOUT_WORLD));
214 PetscCall(DMDestroy(&dm));
215 PetscFunctionReturn(PETSC_SUCCESS);
216 }
217
TestUniversalLabel(MPI_Comm comm)218 static PetscErrorCode TestUniversalLabel(MPI_Comm comm)
219 {
220 DM dm1, dm2;
221 DMLabel bd1, bd2, ulabel;
222 DMUniversalLabel universal;
223 PetscInt pStart, pEnd, p;
224 PetscBool run = PETSC_FALSE, notFile;
225
226 PetscFunctionBeginUser;
227 PetscCall(PetscOptionsGetBool(NULL, NULL, "-universal", &run, NULL));
228 if (!run) PetscFunctionReturn(PETSC_SUCCESS);
229
230 char filename[PETSC_MAX_PATH_LEN];
231 PetscBool flg;
232
233 PetscCall(PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), &flg));
234 if (flg) {
235 PetscCall(DMPlexCreateFromFile(comm, filename, "ex11_plex", PETSC_TRUE, &dm1));
236 } else {
237 PetscCall(DMCreate(comm, &dm1));
238 PetscCall(DMSetType(dm1, DMPLEX));
239 PetscCall(DMSetFromOptions(dm1));
240 }
241 PetscCall(DMHasLabel(dm1, "marker", ¬File));
242 if (notFile) {
243 PetscCall(DMCreateLabel(dm1, "Boundary Faces"));
244 PetscCall(DMGetLabel(dm1, "Boundary Faces", &bd1));
245 PetscCall(DMPlexMarkBoundaryFaces(dm1, 13, bd1));
246 PetscCall(DMCreateLabel(dm1, "Boundary"));
247 PetscCall(DMGetLabel(dm1, "Boundary", &bd2));
248 PetscCall(DMPlexMarkBoundaryFaces(dm1, 121, bd2));
249 PetscCall(DMPlexLabelComplete(dm1, bd2));
250 }
251 PetscCall(PetscObjectSetName((PetscObject)dm1, "First Mesh"));
252 PetscCall(DMViewFromOptions(dm1, NULL, "-dm_view"));
253
254 PetscCall(DMUniversalLabelCreate(dm1, &universal));
255 PetscCall(DMUniversalLabelGetLabel(universal, &ulabel));
256 PetscCall(PetscObjectViewFromOptions((PetscObject)ulabel, NULL, "-universal_view"));
257
258 if (!notFile) {
259 PetscInt Nl, l;
260
261 PetscCall(DMClone(dm1, &dm2));
262 PetscCall(DMGetNumLabels(dm2, &Nl));
263 for (l = Nl - 1; l >= 0; --l) {
264 PetscBool isdepth, iscelltype;
265 const char *name;
266
267 PetscCall(DMGetLabelName(dm2, l, &name));
268 PetscCall(PetscStrncmp(name, "depth", 6, &isdepth));
269 PetscCall(PetscStrncmp(name, "celltype", 9, &iscelltype));
270 if (!isdepth && !iscelltype) PetscCall(DMRemoveLabel(dm2, name, NULL));
271 }
272 } else {
273 PetscCall(DMCreate(comm, &dm2));
274 PetscCall(DMSetType(dm2, DMPLEX));
275 PetscCall(DMSetFromOptions(dm2));
276 }
277 PetscCall(PetscObjectSetName((PetscObject)dm2, "Second Mesh"));
278 PetscCall(DMUniversalLabelCreateLabels(universal, PETSC_TRUE, dm2));
279 PetscCall(DMPlexGetChart(dm2, &pStart, &pEnd));
280 for (p = pStart; p < pEnd; ++p) {
281 PetscInt val;
282
283 PetscCall(DMLabelGetValue(ulabel, p, &val));
284 if (val < 0) continue;
285 PetscCall(DMUniversalLabelSetLabelValue(universal, dm2, PETSC_TRUE, p, val));
286 }
287 PetscCall(DMViewFromOptions(dm2, NULL, "-dm_view"));
288
289 PetscCall(DMUniversalLabelDestroy(&universal));
290 PetscCall(DMDestroy(&dm1));
291 PetscCall(DMDestroy(&dm2));
292 PetscFunctionReturn(PETSC_SUCCESS);
293 }
294
main(int argc,char ** argv)295 int main(int argc, char **argv)
296 {
297 PetscFunctionBeginUser;
298 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
299 /*PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));*/
300 PetscCall(TestInsertion());
301 PetscCall(TestEmptyStrata(PETSC_COMM_WORLD));
302 PetscCall(TestDistribution(PETSC_COMM_WORLD));
303 PetscCall(TestUniversalLabel(PETSC_COMM_WORLD));
304 PetscCall(PetscFinalize());
305 return 0;
306 }
307
308 /*TEST
309
310 test:
311 suffix: 0
312 requires: triangle
313 test:
314 suffix: 1
315 requires: triangle
316 nsize: 2
317 args: -petscpartitioner_type simple
318
319 testset:
320 suffix: gmsh
321 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.msh -petscpartitioner_type simple
322 test:
323 suffix: 1
324 nsize: 1
325 test:
326 suffix: 2
327 nsize: 2
328
329 testset:
330 suffix: exodusii
331 requires: exodusii
332 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/2Dgrd.exo -petscpartitioner_type simple
333 test:
334 suffix: 1
335 nsize: 1
336 test:
337 suffix: 2
338 nsize: 2
339
340 test:
341 suffix: univ
342 requires: triangle
343 args: -universal -dm_view -universal_view
344
345 test:
346 # Note that the labels differ because we have multiply-marked some points during EGADS creation
347 suffix: univ_egads_sphere
348 requires: egads datafilespath
349 args: -universal -dm_plex_filename ${DATAFILESPATH}/meshes/cad/sphere_example.egadslite -dm_view -universal_view
350
351 test:
352 # Note that the labels differ because we have multiply-marked some points during EGADS creation
353 suffix: univ_egads_ball
354 requires: egads ctetgen datafilespath
355 args: -universal -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/sphere_example.egadslite -dm_view -universal_view
356
357 TEST*/
358