Actual source code: ex1.c

  1: static char help[] = "Tests DMLabel operations.\n\n";

  3: #include <petscdm.h>
  4: #include <petscdmplex.h>
  5: #include <petscdmplextransform.h>

  7: PetscErrorCode ViewLabels(DM dm, PetscViewer viewer)
  8: {
  9:   DMLabel     label;
 10:   const char *labelName, *typeName;
 11:   PetscInt    numLabels, l;

 13:   PetscFunctionBegin;
 14:   /* query the number and name of labels*/
 15:   PetscCall(DMGetNumLabels(dm, &numLabels));
 16:   PetscCall(PetscViewerASCIIPrintf(viewer, "Number of labels: %" PetscInt_FMT "\n", numLabels));
 17:   for (l = 0; l < numLabels; ++l) {
 18:     IS labelIS, tmpIS;

 20:     PetscCall(DMGetLabelName(dm, l, &labelName));
 21:     PetscCall(DMGetLabel(dm, labelName, &label));
 22:     PetscCall(DMLabelGetType(label, &typeName));
 23:     PetscCall(PetscViewerASCIIPrintf(viewer, "Label %" PetscInt_FMT ": name: %s type: %s\n", l, labelName, typeName));
 24:     PetscCall(PetscViewerASCIIPrintf(viewer, "IS of values\n"));
 25:     PetscCall(DMLabelGetValueIS(label, &labelIS));
 26:     PetscCall(ISOnComm(labelIS, PetscObjectComm((PetscObject)viewer), PETSC_USE_POINTER, &tmpIS));
 27:     PetscCall(PetscViewerASCIIPushTab(viewer));
 28:     PetscCall(ISView(tmpIS, viewer));
 29:     PetscCall(PetscViewerASCIIPopTab(viewer));
 30:     PetscCall(ISDestroy(&tmpIS));
 31:     PetscCall(ISDestroy(&labelIS));
 32:     PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
 33:   }
 34:   /* Making sure that string literals work */
 35:   PetscCall(PetscViewerASCIIPrintf(viewer, "\n\nCell Set label IS\n"));
 36:   PetscCall(DMGetLabel(dm, "Cell Sets", &label));
 37:   if (label) {
 38:     IS labelIS, tmpIS;

 40:     PetscCall(DMLabelGetValueIS(label, &labelIS));
 41:     PetscCall(ISOnComm(labelIS, PetscObjectComm((PetscObject)viewer), PETSC_USE_POINTER, &tmpIS));
 42:     PetscCall(ISView(tmpIS, viewer));
 43:     PetscCall(ISDestroy(&tmpIS));
 44:     PetscCall(ISDestroy(&labelIS));
 45:   }
 46:   PetscFunctionReturn(PETSC_SUCCESS);
 47: }

 49: PetscErrorCode CheckLabelsSame(DMLabel label0, DMLabel label1)
 50: {
 51:   const char *name0, *name1;
 52:   PetscBool   same;
 53:   char       *msg;

 55:   PetscFunctionBegin;
 56:   PetscCall(PetscObjectGetName((PetscObject)label0, &name0));
 57:   PetscCall(PetscObjectGetName((PetscObject)label1, &name1));
 58:   PetscCall(DMLabelCompare(PETSC_COMM_WORLD, label0, label1, &same, &msg));
 59:   PetscCheck(same == (PetscBool)!msg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "DMLabelCompare returns inconsistent same=%d msg=\"%s\"", same, msg);
 60:   PetscCheck(same, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Labels \"%s\" and \"%s\" should not differ! Message:\n%s", name0, name1, msg);
 61:   /* Test passing NULL, must not fail */
 62:   PetscCall(DMLabelCompare(PETSC_COMM_WORLD, label0, label1, NULL, NULL));
 63:   PetscCall(PetscFree(msg));
 64:   PetscFunctionReturn(PETSC_SUCCESS);
 65: }

 67: PetscErrorCode CheckLabelsNotSame(DMLabel label0, DMLabel label1)
 68: {
 69:   const char *name0, *name1;
 70:   PetscBool   same;
 71:   char       *msg;

 73:   PetscFunctionBegin;
 74:   PetscCall(PetscObjectGetName((PetscObject)label0, &name0));
 75:   PetscCall(PetscObjectGetName((PetscObject)label1, &name1));
 76:   PetscCall(DMLabelCompare(PETSC_COMM_WORLD, label0, label1, &same, &msg));
 77:   PetscCheck(same == (PetscBool)!msg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "DMLabelCompare returns inconsistent same=%d msg=\"%s\"", same, msg);
 78:   PetscCheck(!same, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Labels \"%s\" and \"%s\" should differ!", name0, name1);
 79:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Compare label \"%s\" with \"%s\": %s\n", name0, name1, msg));
 80:   PetscCall(PetscFree(msg));
 81:   PetscFunctionReturn(PETSC_SUCCESS);
 82: }

 84: PetscErrorCode CheckDMLabelsSame(DM dm0, DM dm1)
 85: {
 86:   const char *name0, *name1;
 87:   PetscBool   same;
 88:   char       *msg;

 90:   PetscFunctionBegin;
 91:   PetscCall(PetscObjectGetName((PetscObject)dm0, &name0));
 92:   PetscCall(PetscObjectGetName((PetscObject)dm1, &name1));
 93:   PetscCall(DMCompareLabels(dm0, dm1, &same, &msg));
 94:   PetscCheck(same == (PetscBool)!msg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "DMCompareLabels returns inconsistent same=%d msg=\"%s\"", same, msg);
 95:   PetscCheck(same, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Labels of DMs \"%s\" and \"%s\" should not differ! Message:\n%s", name0, name1, msg);
 96:   /* Test passing NULL, must not fail */
 97:   PetscCall(DMCompareLabels(dm0, dm1, NULL, NULL));
 98:   PetscCall(PetscFree(msg));
 99:   PetscFunctionReturn(PETSC_SUCCESS);
100: }

102: PetscErrorCode CheckDMLabelsNotSame(DM dm0, DM dm1)
103: {
104:   const char *name0, *name1;
105:   PetscBool   same;
106:   char       *msg;

108:   PetscFunctionBegin;
109:   PetscCall(PetscObjectGetName((PetscObject)dm0, &name0));
110:   PetscCall(PetscObjectGetName((PetscObject)dm1, &name1));
111:   PetscCall(DMCompareLabels(dm0, dm1, &same, &msg));
112:   PetscCheck(same == (PetscBool)!msg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "DMCompareLabels returns inconsistent same=%d msg=\"%s\"", same, msg);
113:   PetscCheck(!same, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Labels of DMs \"%s\" and \"%s\" should differ!", name0, name1);
114:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Labels of DMs \"%s\" and \"%s\" differ: %s\n", name0, name1, msg));
115:   PetscCall(PetscFree(msg));
116:   PetscFunctionReturn(PETSC_SUCCESS);
117: }

119: PetscErrorCode CreateMesh(const char name[], DM *newdm)
120: {
121:   DM        dm, dmDist;
122:   char      filename[PETSC_MAX_PATH_LEN] = "";
123:   PetscBool interpolate                  = PETSC_FALSE;

125:   PetscFunctionBegin;
126:   /* initialize and get options */
127:   PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "DMLabel ex1 Options", "DMLabel");
128:   PetscCall(PetscOptionsString("-i", "filename to read", "ex1.c", filename, filename, sizeof(filename), NULL));
129:   PetscCall(PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex1.c", interpolate, &interpolate, NULL));
130:   PetscOptionsEnd();

132:   /* create and distribute DM */
133:   PetscCall(DMPlexCreateFromFile(PETSC_COMM_WORLD, filename, "ex1_plex", interpolate, &dm));
134:   PetscCall(DMPlexDistribute(dm, 0, NULL, &dmDist));
135:   if (dmDist) {
136:     PetscCall(DMDestroy(&dm));
137:     dm = dmDist;
138:   }
139:   PetscCall(DMSetFromOptions(dm));
140:   PetscCall(PetscObjectSetName((PetscObject)dm, name));
141:   *newdm = dm;
142:   PetscFunctionReturn(PETSC_SUCCESS);
143: }

145: static PetscErrorCode TestEphemeralLabels(DM dm)
146: {
147:   DMPlexTransform tr;
148:   DM              tdm;
149:   DMLabel         label, labelTmp;

151:   PetscFunctionBeginUser;
152:   PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), &tr));
153:   PetscCall(PetscObjectSetName((PetscObject)tr, "Transform"));
154:   PetscCall(DMPlexTransformSetDM(tr, dm));
155:   PetscCall(DMPlexTransformSetFromOptions(tr));
156:   PetscCall(DMPlexTransformSetUp(tr));

158:   PetscCall(DMPlexCreateEphemeral(tr, &tdm));
159:   PetscCall(DMPlexTransformDestroy(&tr));
160:   PetscCall(PetscObjectSetName((PetscObject)tdm, "Ephemeral Mesh"));
161:   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)tdm, "eph_"));

163:   PetscCall(DMGetLabel(tdm, "marker", &label));
164:   PetscCall(DMLabelDuplicate(label, &labelTmp));
165:   PetscCall(CheckLabelsSame(label, labelTmp));
166:   PetscCall(DMLabelDestroy(&labelTmp));
167:   PetscCall(DMDestroy(&tdm));
168:   PetscFunctionReturn(PETSC_SUCCESS);
169: }

171: int main(int argc, char **argv)
172: {
173:   DM dm;

175:   PetscFunctionBeginUser;
176:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
177:   PetscCall(CreateMesh("plex0", &dm));
178:   /* add custom labels to test adding/removal */
179:   {
180:     DMLabel  label0, label1, label2, label3;
181:     PetscInt p, pStart, pEnd;
182:     PetscCall(DMPlexGetChart(dm, &pStart, &pEnd));
183:     /* create label in DM and get from DM */
184:     PetscCall(DMCreateLabel(dm, "label0"));
185:     PetscCall(DMGetLabel(dm, "label0", &label0));
186:     /* alternative: create standalone label and add to DM; needs to be destroyed */
187:     PetscCall(DMLabelCreate(PETSC_COMM_SELF, "label1", &label1));
188:     PetscCall(DMAddLabel(dm, label1));

190:     pEnd = PetscMin(pEnd, pStart + 5);
191:     for (p = pStart; p < pEnd; p++) {
192:       PetscCall(DMLabelSetValue(label0, p, 1));
193:       PetscCall(DMLabelSetValue(label1, p, 2));
194:     }
195:     /* duplicate label */
196:     PetscCall(DMLabelDuplicate(label0, &label2));
197:     PetscCall(DMLabelDuplicate(label1, &label3));
198:     PetscCall(PetscObjectSetName((PetscObject)label2, "label2"));
199:     PetscCall(PetscObjectSetName((PetscObject)label3, "label3"));
200:     PetscCall(DMAddLabel(dm, label2));
201:     PetscCall(DMAddLabel(dm, label3));
202:     /* remove the labels in this scope */
203:     PetscCall(DMLabelDestroy(&label1));
204:     PetscCall(DMLabelDestroy(&label2));
205:     PetscCall(DMLabelDestroy(&label3));
206:   }

208:   PetscCall(ViewLabels(dm, PETSC_VIEWER_STDOUT_WORLD));

210:   /* do label perturbations and comparisons */
211:   {
212:     DMLabel  label0, label1, label2, label3;
213:     PetscInt val;
214:     PetscInt p, pStart, pEnd;

216:     PetscCall(DMGetLabel(dm, "label0", &label0));
217:     PetscCall(DMGetLabel(dm, "label1", &label1));
218:     PetscCall(DMGetLabel(dm, "label2", &label2));
219:     PetscCall(DMGetLabel(dm, "label3", &label3));

221:     PetscCall(CheckLabelsNotSame(label0, label1));
222:     PetscCall(CheckLabelsSame(label0, label2));
223:     PetscCall(CheckLabelsSame(label1, label3));

225:     PetscCall(DMLabelGetDefaultValue(label1, &val));
226:     PetscCall(DMLabelSetDefaultValue(label1, 333));
227:     PetscCall(CheckLabelsNotSame(label1, label3));
228:     PetscCall(DMLabelSetDefaultValue(label1, val));
229:     PetscCall(CheckLabelsSame(label1, label3));

231:     PetscCall(DMLabelGetBounds(label1, &pStart, &pEnd));

233:     for (p = pStart; p < pEnd; p++) {
234:       PetscCall(DMLabelGetValue(label1, p, &val));
235:       // This is weird. Perhaps we should not need to call DMLabelClearValue()
236:       PetscCall(DMLabelClearValue(label1, p, val));
237:       val++;
238:       PetscCall(DMLabelSetValue(label1, p, val));
239:     }
240:     PetscCall(CheckLabelsNotSame(label1, label3));
241:     for (p = pStart; p < pEnd; p++) {
242:       PetscCall(DMLabelGetValue(label1, p, &val));
243:       // This is weird. Perhaps we should not need to call DMLabelClearValue()
244:       PetscCall(DMLabelClearValue(label1, p, val));
245:       val--;
246:       PetscCall(DMLabelSetValue(label1, p, val));
247:     }
248:     PetscCall(CheckLabelsSame(label1, label3));

250:     PetscCall(DMLabelGetValue(label3, pEnd - 1, &val));
251:     PetscCall(DMLabelSetValue(label3, pEnd, val));
252:     PetscCall(CheckLabelsNotSame(label1, label3));
253:     // This is weird. Perhaps we should not need to call DMLabelClearValue()
254:     PetscCall(DMLabelClearValue(label3, pEnd, val));
255:     PetscCall(CheckLabelsSame(label1, label3));
256:   }

258:   {
259:     DM       dm1;
260:     DMLabel  label02, label12;
261:     PetscInt p = 0, val;

263:     PetscCall(CreateMesh("plex1", &dm1));
264:     PetscCall(CheckDMLabelsNotSame(dm, dm1));

266:     PetscCall(DMCopyLabels(dm, dm1, PETSC_OWN_POINTER, PETSC_FALSE, DM_COPY_LABELS_REPLACE));
267:     PetscCall(CheckDMLabelsSame(dm, dm1));

269:     PetscCall(DMCopyLabels(dm, dm1, PETSC_COPY_VALUES, PETSC_FALSE, DM_COPY_LABELS_REPLACE));
270:     PetscCall(DMGetLabel(dm, "label2", &label02));
271:     PetscCall(DMGetLabel(dm1, "label2", &label12));
272:     PetscCall(CheckLabelsSame(label02, label12));

274:     PetscCall(DMLabelGetValue(label12, p, &val));
275:     // This is weird. Perhaps we should not need to call DMLabelClearValue()
276:     PetscCall(DMLabelClearValue(label12, p, val));
277:     PetscCall(DMLabelSetValue(label12, p, val + 1));
278:     PetscCall(CheckLabelsNotSame(label02, label12));
279:     PetscCall(CheckDMLabelsNotSame(dm, dm1));

281:     // This is weird. Perhaps we should not need to call DMLabelClearValue()
282:     PetscCall(DMLabelClearValue(label12, p, val + 1));
283:     PetscCall(DMLabelSetValue(label12, p, val));
284:     PetscCall(CheckLabelsSame(label02, label12));
285:     PetscCall(CheckDMLabelsSame(dm, dm1));

287:     PetscCall(PetscObjectSetName((PetscObject)label12, "label12"));
288:     PetscCall(CheckDMLabelsNotSame(dm, dm1));
289:     PetscCall(PetscObjectSetName((PetscObject)label12, "label2"));
290:     PetscCall(CheckDMLabelsSame(dm, dm1));

292:     PetscCall(DMDestroy(&dm1));
293:   }
294:   // Test adding strata and filtering
295:   {
296:     DMLabel  labelA, labelB;
297:     IS       valueIS;
298:     PetscInt pStart, pEnd, lStart = 5, lEnd = 19;

300:     PetscCall(DMPlexGetChart(dm, &pStart, &pEnd));
301:     PetscCall(DMCreateLabel(dm, "labelA"));
302:     PetscCall(DMCreateLabel(dm, "labelB"));
303:     PetscCall(DMGetLabel(dm, "labelA", &labelA));
304:     PetscCall(DMGetLabel(dm, "labelB", &labelB));
305:     for (PetscInt p = pStart; p < pEnd; ++p) {
306:       if (p < lStart || p >= lEnd) continue;
307:       if (p % 2) PetscCall(DMLabelSetValue(labelA, p, 19));
308:       else PetscCall(DMLabelSetValue(labelA, p, 17));
309:     }
310:     PetscCall(DMLabelGetValueIS(labelA, &valueIS));
311:     PetscCall(DMLabelAddStrataIS(labelA, valueIS));
312:     PetscCall(ISDestroy(&valueIS));
313:     for (PetscInt p = pStart; p < pEnd; ++p) {
314:       if (p % 2) PetscCall(DMLabelSetValue(labelB, p, 19));
315:       else PetscCall(DMLabelSetValue(labelB, p, 17));
316:     }
317:     PetscCall(DMLabelFilter(labelB, lStart, lEnd));
318:     PetscCall(CheckLabelsSame(labelA, labelB));
319:     PetscCall(DMRemoveLabel(dm, "labelA", NULL));
320:     PetscCall(DMRemoveLabel(dm, "labelB", NULL));
321:   }

323:   /* remove label0 and label1 just to test manual removal; let label3 be removed automatically by DMDestroy() */
324:   {
325:     DMLabel label0, label1, label2;
326:     PetscCall(DMGetLabel(dm, "label0", &label0));
327:     PetscCall(DMGetLabel(dm, "label1", &label1));
328:     PetscCheck(label0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label0 must not be NULL now");
329:     PetscCheck(label1, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label1 must not be NULL now");
330:     PetscCall(DMRemoveLabel(dm, "label1", NULL));
331:     PetscCall(DMRemoveLabel(dm, "label2", &label2));
332:     PetscCall(DMRemoveLabelBySelf(dm, &label0, PETSC_TRUE));
333:     PetscCall(DMGetLabel(dm, "label0", &label0));
334:     PetscCall(DMGetLabel(dm, "label1", &label1));
335:     PetscCheck(!label0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label0 must be NULL now");
336:     PetscCheck(!label1, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label1 must be NULL now");
337:     PetscCheck(label2, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label2 must not be NULL now");
338:     PetscCall(DMRemoveLabelBySelf(dm, &label2, PETSC_FALSE)); /* this should do nothing */
339:     PetscCheck(label2, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label2 must not be NULL now");
340:     PetscCall(DMLabelDestroy(&label2));
341:     PetscCall(DMGetLabel(dm, "label2", &label2));
342:     PetscCheck(!label2, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "label2 must be NULL now");
343:   }

345:   PetscCall(TestEphemeralLabels(dm));

347:   PetscCall(DMDestroy(&dm));
348:   PetscCall(PetscFinalize());
349:   return 0;
350: }

352: /*TEST

354:   test:
355:     suffix: 0
356:     nsize: {{1 2}separate output}
357:     args: -i ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo -interpolate
358:     requires: exodusii

360: TEST*/