Actual source code: ex97.c
1: static char help[] = "Test DMPlexGetCellType\n\n";
3: #include <petsc.h>
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;
14: PetscInitialize(&argc, &argv, NULL, help);
15: PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "FEM Layout Options", "ex97");
16: PetscOptionsString("-i", "Filename to read", "ex97", ifilename, ifilename, sizeof(ifilename), NULL);
17: PetscOptionsEnd();
19: DMPlexCreateFromFile(PETSC_COMM_WORLD, ifilename, NULL, PETSC_TRUE, &dm);
20: DMPlexDistributeSetDefault(dm, PETSC_FALSE);
21: DMSetFromOptions(dm);
23: DMPlexDistribute(dm, 0, NULL, &pdm);
24: if (pdm) {
25: DMDestroy(&dm);
26: dm = pdm;
27: }
28: PetscObjectSetName((PetscObject)dm, "ex97");
29: DMViewFromOptions(dm, NULL, "-dm_view");
31: DMGetLabel(dm, "celltype", &label);
32: DMLabelView(label, PETSC_VIEWER_STDOUT_WORLD);
33: DMPlexGetHeightStratum(dm, 0, &pStart, &pEnd);
34: for (p = pStart; p < pEnd; ++p) {
35: DMPlexGetCellType(dm, p, &cellType);
36: PetscPrintf(PETSC_COMM_SELF, "cell: %" PetscInt_FMT " type: %d\n", p, cellType);
37: }
38: DMDestroy(&dm);
40: PetscFinalize();
41: return 0;
42: }
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*/