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: %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: %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, "eph_", &tdm));
    159:   PetscCall(DMPlexTransformDestroy(&tr));
    160:   PetscCall(PetscObjectSetName((PetscObject)tdm, "Ephemeral Mesh"));

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

    344:   PetscCall(TestEphemeralLabels(dm));

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

    351: /*TEST

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

    359: TEST*/