xref: /petsc/src/dm/label/tutorials/ex1.c (revision 58d68138c660dfb4e9f5b03334792cd4f2ffd7cc)
1 static char help[] = "Tests DMLabel operations.\n\n";
2 
3 #include <petscdm.h>
4 #include <petscdmplex.h>
5 
6 PetscErrorCode ViewLabels(DM dm, PetscViewer viewer) {
7   DMLabel     label;
8   const char *labelName;
9   PetscInt    numLabels, l;
10 
11   PetscFunctionBegin;
12   /* query the number and name of labels*/
13   PetscCall(DMGetNumLabels(dm, &numLabels));
14   PetscCall(PetscViewerASCIIPrintf(viewer, "Number of labels: %" PetscInt_FMT "\n", numLabels));
15   for (l = 0; l < numLabels; ++l) {
16     IS labelIS, tmpIS;
17 
18     PetscCall(DMGetLabelName(dm, l, &labelName));
19     PetscCall(PetscViewerASCIIPrintf(viewer, "Label %" PetscInt_FMT ": name: %s\n", l, labelName));
20     PetscCall(PetscViewerASCIIPrintf(viewer, "IS of values\n"));
21     PetscCall(DMGetLabel(dm, labelName, &label));
22     PetscCall(DMLabelGetValueIS(label, &labelIS));
23     PetscCall(ISOnComm(labelIS, PetscObjectComm((PetscObject)viewer), PETSC_USE_POINTER, &tmpIS));
24     PetscCall(PetscViewerASCIIPushTab(viewer));
25     PetscCall(ISView(tmpIS, viewer));
26     PetscCall(PetscViewerASCIIPopTab(viewer));
27     PetscCall(ISDestroy(&tmpIS));
28     PetscCall(ISDestroy(&labelIS));
29     PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
30   }
31   /* Making sure that string literals work */
32   PetscCall(PetscViewerASCIIPrintf(viewer, "\n\nCell Set label IS\n"));
33   PetscCall(DMGetLabel(dm, "Cell Sets", &label));
34   if (label) {
35     IS labelIS, tmpIS;
36 
37     PetscCall(DMLabelGetValueIS(label, &labelIS));
38     PetscCall(ISOnComm(labelIS, PetscObjectComm((PetscObject)viewer), PETSC_USE_POINTER, &tmpIS));
39     PetscCall(ISView(tmpIS, viewer));
40     PetscCall(ISDestroy(&tmpIS));
41     PetscCall(ISDestroy(&labelIS));
42   }
43   PetscFunctionReturn(0);
44 }
45 
46 PetscErrorCode CheckLabelsSame(DMLabel label0, DMLabel label1) {
47   const char *name0, *name1;
48   PetscBool   same;
49   char       *msg;
50 
51   PetscFunctionBegin;
52   PetscCall(PetscObjectGetName((PetscObject)label0, &name0));
53   PetscCall(PetscObjectGetName((PetscObject)label1, &name1));
54   PetscCall(DMLabelCompare(PETSC_COMM_WORLD, label0, label1, &same, &msg));
55   PetscCheck(same == (PetscBool)!msg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "DMLabelCompare returns inconsistent same=%d msg=\"%s\"", same, msg);
56   PetscCheck(same, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Labels \"%s\" and \"%s\" should not differ! Message:\n%s", name0, name1, msg);
57   /* Test passing NULL, must not fail */
58   PetscCall(DMLabelCompare(PETSC_COMM_WORLD, label0, label1, NULL, NULL));
59   PetscCall(PetscFree(msg));
60   PetscFunctionReturn(0);
61 }
62 
63 PetscErrorCode CheckLabelsNotSame(DMLabel label0, DMLabel label1) {
64   const char *name0, *name1;
65   PetscBool   same;
66   char       *msg;
67 
68   PetscFunctionBegin;
69   PetscCall(PetscObjectGetName((PetscObject)label0, &name0));
70   PetscCall(PetscObjectGetName((PetscObject)label1, &name1));
71   PetscCall(DMLabelCompare(PETSC_COMM_WORLD, label0, label1, &same, &msg));
72   PetscCheck(same == (PetscBool)!msg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "DMLabelCompare returns inconsistent same=%d msg=\"%s\"", same, msg);
73   PetscCheck(!same, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Labels \"%s\" and \"%s\" should differ!", name0, name1);
74   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Compare label \"%s\" with \"%s\": %s\n", name0, name1, msg));
75   PetscCall(PetscFree(msg));
76   PetscFunctionReturn(0);
77 }
78 
79 PetscErrorCode CheckDMLabelsSame(DM dm0, DM dm1) {
80   const char *name0, *name1;
81   PetscBool   same;
82   char       *msg;
83 
84   PetscFunctionBegin;
85   PetscCall(PetscObjectGetName((PetscObject)dm0, &name0));
86   PetscCall(PetscObjectGetName((PetscObject)dm1, &name1));
87   PetscCall(DMCompareLabels(dm0, dm1, &same, &msg));
88   PetscCheck(same == (PetscBool)!msg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "DMCompareLabels returns inconsistent same=%d msg=\"%s\"", same, msg);
89   PetscCheck(same, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Labels of DMs \"%s\" and \"%s\" should not differ! Message:\n%s", name0, name1, msg);
90   /* Test passing NULL, must not fail */
91   PetscCall(DMCompareLabels(dm0, dm1, NULL, NULL));
92   PetscCall(PetscFree(msg));
93   PetscFunctionReturn(0);
94 }
95 
96 PetscErrorCode CheckDMLabelsNotSame(DM dm0, DM dm1) {
97   const char *name0, *name1;
98   PetscBool   same;
99   char       *msg;
100 
101   PetscFunctionBegin;
102   PetscCall(PetscObjectGetName((PetscObject)dm0, &name0));
103   PetscCall(PetscObjectGetName((PetscObject)dm1, &name1));
104   PetscCall(DMCompareLabels(dm0, dm1, &same, &msg));
105   PetscCheck(same == (PetscBool)!msg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "DMCompareLabels returns inconsistent same=%d msg=\"%s\"", same, msg);
106   PetscCheck(!same, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Labels of DMs \"%s\" and \"%s\" should differ!", name0, name1);
107   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Labels of DMs \"%s\" and \"%s\" differ: %s\n", name0, name1, msg));
108   PetscCall(PetscFree(msg));
109   PetscFunctionReturn(0);
110 }
111 
112 PetscErrorCode CreateMesh(const char name[], DM *newdm) {
113   DM        dm, dmDist;
114   char      filename[PETSC_MAX_PATH_LEN] = "";
115   PetscBool interpolate                  = PETSC_FALSE;
116 
117   PetscFunctionBegin;
118   /* initialize and get options */
119   PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "DMLabel ex1 Options", "DMLabel");
120   PetscCall(PetscOptionsString("-i", "filename to read", "ex1.c", filename, filename, sizeof(filename), NULL));
121   PetscCall(PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex1.c", interpolate, &interpolate, NULL));
122   PetscOptionsEnd();
123 
124   /* create and distribute DM */
125   PetscCall(DMPlexCreateFromFile(PETSC_COMM_WORLD, filename, "ex1_plex", interpolate, &dm));
126   PetscCall(DMPlexDistribute(dm, 0, NULL, &dmDist));
127   if (dmDist) {
128     PetscCall(DMDestroy(&dm));
129     dm = dmDist;
130   }
131   PetscCall(DMSetFromOptions(dm));
132   PetscCall(PetscObjectSetName((PetscObject)dm, name));
133   *newdm = dm;
134   PetscFunctionReturn(0);
135 }
136 
137 int main(int argc, char **argv) {
138   DM dm;
139 
140   PetscFunctionBeginUser;
141   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
142   PetscCall(CreateMesh("plex0", &dm));
143   /* add custom labels to test adding/removal */
144   {
145     DMLabel  label0, label1, label2, label3;
146     PetscInt p, pStart, pEnd;
147     PetscCall(DMPlexGetChart(dm, &pStart, &pEnd));
148     /* create label in DM and get from DM */
149     PetscCall(DMCreateLabel(dm, "label0"));
150     PetscCall(DMGetLabel(dm, "label0", &label0));
151     /* alternative: create standalone label and add to DM; needs to be destroyed */
152     PetscCall(DMLabelCreate(PETSC_COMM_SELF, "label1", &label1));
153     PetscCall(DMAddLabel(dm, label1));
154 
155     pEnd = PetscMin(pEnd, pStart + 5);
156     for (p = pStart; p < pEnd; p++) {
157       PetscCall(DMLabelSetValue(label0, p, 1));
158       PetscCall(DMLabelSetValue(label1, p, 2));
159     }
160     /* duplicate label */
161     PetscCall(DMLabelDuplicate(label0, &label2));
162     PetscCall(DMLabelDuplicate(label1, &label3));
163     PetscCall(PetscObjectSetName((PetscObject)label2, "label2"));
164     PetscCall(PetscObjectSetName((PetscObject)label3, "label3"));
165     PetscCall(DMAddLabel(dm, label2));
166     PetscCall(DMAddLabel(dm, label3));
167     /* remove the labels in this scope */
168     PetscCall(DMLabelDestroy(&label1));
169     PetscCall(DMLabelDestroy(&label2));
170     PetscCall(DMLabelDestroy(&label3));
171   }
172 
173   PetscCall(ViewLabels(dm, PETSC_VIEWER_STDOUT_WORLD));
174 
175   /* do label perturbations and comparisons */
176   {
177     DMLabel  label0, label1, label2, label3;
178     PetscInt val;
179     PetscInt p, pStart, pEnd;
180 
181     PetscCall(DMGetLabel(dm, "label0", &label0));
182     PetscCall(DMGetLabel(dm, "label1", &label1));
183     PetscCall(DMGetLabel(dm, "label2", &label2));
184     PetscCall(DMGetLabel(dm, "label3", &label3));
185 
186     PetscCall(CheckLabelsNotSame(label0, label1));
187     PetscCall(CheckLabelsSame(label0, label2));
188     PetscCall(CheckLabelsSame(label1, label3));
189 
190     PetscCall(DMLabelGetDefaultValue(label1, &val));
191     PetscCall(DMLabelSetDefaultValue(label1, 333));
192     PetscCall(CheckLabelsNotSame(label1, label3));
193     PetscCall(DMLabelSetDefaultValue(label1, val));
194     PetscCall(CheckLabelsSame(label1, label3));
195 
196     PetscCall(DMLabelGetBounds(label1, &pStart, &pEnd));
197 
198     for (p = pStart; p < pEnd; p++) {
199       PetscCall(DMLabelGetValue(label1, p, &val));
200       // This is weird. Perhaps we should not need to call DMLabelClearValue()
201       PetscCall(DMLabelClearValue(label1, p, val));
202       val++;
203       PetscCall(DMLabelSetValue(label1, p, val));
204     }
205     PetscCall(CheckLabelsNotSame(label1, label3));
206     for (p = pStart; p < pEnd; p++) {
207       PetscCall(DMLabelGetValue(label1, p, &val));
208       // This is weird. Perhaps we should not need to call DMLabelClearValue()
209       PetscCall(DMLabelClearValue(label1, p, val));
210       val--;
211       PetscCall(DMLabelSetValue(label1, p, val));
212     }
213     PetscCall(CheckLabelsSame(label1, label3));
214 
215     PetscCall(DMLabelGetValue(label3, pEnd - 1, &val));
216     PetscCall(DMLabelSetValue(label3, pEnd, val));
217     PetscCall(CheckLabelsNotSame(label1, label3));
218     // This is weird. Perhaps we should not need to call DMLabelClearValue()
219     PetscCall(DMLabelClearValue(label3, pEnd, val));
220     PetscCall(CheckLabelsSame(label1, label3));
221   }
222 
223   {
224     DM       dm1;
225     DMLabel  label02, label12;
226     PetscInt p = 0, val;
227 
228     PetscCall(CreateMesh("plex1", &dm1));
229     PetscCall(CheckDMLabelsNotSame(dm, dm1));
230 
231     PetscCall(DMCopyLabels(dm, dm1, PETSC_OWN_POINTER, PETSC_FALSE, DM_COPY_LABELS_REPLACE));
232     PetscCall(CheckDMLabelsSame(dm, dm1));
233 
234     PetscCall(DMCopyLabels(dm, dm1, PETSC_COPY_VALUES, PETSC_FALSE, DM_COPY_LABELS_REPLACE));
235     PetscCall(DMGetLabel(dm, "label2", &label02));
236     PetscCall(DMGetLabel(dm1, "label2", &label12));
237     PetscCall(CheckLabelsSame(label02, label12));
238 
239     PetscCall(DMLabelGetValue(label12, p, &val));
240     // This is weird. Perhaps we should not need to call DMLabelClearValue()
241     PetscCall(DMLabelClearValue(label12, p, val));
242     PetscCall(DMLabelSetValue(label12, p, val + 1));
243     PetscCall(CheckLabelsNotSame(label02, label12));
244     PetscCall(CheckDMLabelsNotSame(dm, dm1));
245 
246     // This is weird. Perhaps we should not need to call DMLabelClearValue()
247     PetscCall(DMLabelClearValue(label12, p, val + 1));
248     PetscCall(DMLabelSetValue(label12, p, val));
249     PetscCall(CheckLabelsSame(label02, label12));
250     PetscCall(CheckDMLabelsSame(dm, dm1));
251 
252     PetscCall(PetscObjectSetName((PetscObject)label12, "label12"));
253     PetscCall(CheckDMLabelsNotSame(dm, dm1));
254     PetscCall(PetscObjectSetName((PetscObject)label12, "label2"));
255     PetscCall(CheckDMLabelsSame(dm, dm1));
256 
257     PetscCall(DMDestroy(&dm1));
258   }
259 
260   /* remove label0 and label1 just to test manual removal; let label3 be removed automatically by DMDestroy() */
261   {
262     DMLabel label0, label1, label2;
263     PetscCall(DMGetLabel(dm, "label0", &label0));
264     PetscCall(DMGetLabel(dm, "label1", &label1));
265     PetscCheck(label0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label0 must not be NULL now");
266     PetscCheck(label1, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label1 must not be NULL now");
267     PetscCall(DMRemoveLabel(dm, "label1", NULL));
268     PetscCall(DMRemoveLabel(dm, "label2", &label2));
269     PetscCall(DMRemoveLabelBySelf(dm, &label0, PETSC_TRUE));
270     PetscCall(DMGetLabel(dm, "label0", &label0));
271     PetscCall(DMGetLabel(dm, "label1", &label1));
272     PetscCheck(!label0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label0 must be NULL now");
273     PetscCheck(!label1, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label1 must be NULL now");
274     PetscCheck(label2, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label2 must not be NULL now");
275     PetscCall(DMRemoveLabelBySelf(dm, &label2, PETSC_FALSE)); /* this should do nothing */
276     PetscCheck(label2, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label2 must not be NULL now");
277     PetscCall(DMLabelDestroy(&label2));
278     PetscCall(DMGetLabel(dm, "label2", &label2));
279     PetscCheck(!label2, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label2 must be NULL now");
280   }
281 
282   PetscCall(DMDestroy(&dm));
283   PetscCall(PetscFinalize());
284   return 0;
285 }
286 
287 /*TEST
288 
289   test:
290     suffix: 0
291     nsize: {{1 2}separate output}
292     args: -i ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo -interpolate
293     requires: exodusii
294 
295 TEST*/
296