1 static char help[] = "Test DMPlexGetCellType\n\n"; 2 3 #include <petsc.h> 4 5 int main(int argc, char **argv) { 6 DM dm, pdm; 7 char ifilename[PETSC_MAX_PATH_LEN]; 8 PetscInt pStart, pEnd, p; 9 DMPolytopeType cellType; 10 DMLabel label; 11 12 PetscFunctionBeginUser; 13 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 14 PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "FEM Layout Options", "ex97"); 15 PetscCall(PetscOptionsString("-i", "Filename to read", "ex97", ifilename, ifilename, sizeof(ifilename), NULL)); 16 PetscOptionsEnd(); 17 18 PetscCall(DMPlexCreateFromFile(PETSC_COMM_WORLD, ifilename, NULL, PETSC_TRUE, &dm)); 19 PetscCall(DMPlexDistributeSetDefault(dm, PETSC_FALSE)); 20 PetscCall(DMSetFromOptions(dm)); 21 22 PetscCall(DMPlexDistribute(dm, 0, NULL, &pdm)); 23 if (pdm) { 24 PetscCall(DMDestroy(&dm)); 25 dm = pdm; 26 } 27 PetscCall(PetscObjectSetName((PetscObject)dm, "ex97")); 28 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view")); 29 30 PetscCall(DMGetLabel(dm, "celltype", &label)); 31 PetscCall(DMLabelView(label, PETSC_VIEWER_STDOUT_WORLD)); 32 PetscCall(DMPlexGetHeightStratum(dm, 0, &pStart, &pEnd)); 33 for (p = pStart; p < pEnd; ++p) { 34 PetscCall(DMPlexGetCellType(dm, p, &cellType)); 35 PetscCall(PetscPrintf(PETSC_COMM_SELF, "cell: %" PetscInt_FMT " type: %d\n", p, cellType)); 36 } 37 PetscCall(DMDestroy(&dm)); 38 39 PetscCall(PetscFinalize()); 40 return 0; 41 } 42 43 /*TEST 44 build: 45 requires: !complex 46 testset: 47 args: -i ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -dm_view 48 nsize: 1 49 test: 50 suffix: 0 51 args: 52 TEST*/ 53