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