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