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 TestInsertion()
7: {
8: DMLabel label, label2;
9: const PetscInt values[5] = {0, 3, 4, -1, 176}, N = 10000;
10: PetscInt i, v;
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;
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;
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;
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;
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: }
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;
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;
121: PetscCall(DMCreateLabel(dm, "depth"));
122: PetscCall(DMPlexGetDepthLabel(dm, &label));
123: if (rank == 0) {
124: PetscInt i;
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;
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: }
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;
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: }
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;
226: PetscFunctionBeginUser;
227: PetscCall(PetscOptionsGetBool(NULL, NULL, "-universal", &run, NULL));
228: if (!run) PetscFunctionReturn(PETSC_SUCCESS);
230: char filename[PETSC_MAX_PATH_LEN];
231: PetscBool flg;
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"));
254: PetscCall(DMUniversalLabelCreate(dm1, &universal));
255: PetscCall(DMUniversalLabelGetLabel(universal, &ulabel));
256: PetscCall(PetscObjectViewFromOptions((PetscObject)ulabel, NULL, "-universal_view"));
258: if (!notFile) {
259: PetscInt Nl, l;
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;
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;
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"));
289: PetscCall(DMUniversalLabelDestroy(&universal));
290: PetscCall(DMDestroy(&dm1));
291: PetscCall(DMDestroy(&dm2));
292: PetscFunctionReturn(PETSC_SUCCESS);
293: }
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: }
308: /*TEST
310: test:
311: suffix: 0
312: requires: triangle
313: test:
314: suffix: 1
315: requires: triangle
316: nsize: 2
317: args: -petscpartitioner_type simple
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
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
340: test:
341: suffix: univ
342: requires: triangle
343: args: -universal -dm_view -universal_view
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
349: args: -universal -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.egadslite -dm_view -universal_view
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
355: args: -universal -dm_plex_boundary_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.egadslite -dm_view -universal_view
357: TEST*/