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