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