1c4762a1bSJed Brown static char help[] = "Element closure restrictions in tensor/lexicographic/spectral-element ordering using DMPlex\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscdmplex.h> 4c4762a1bSJed Brown 5c4762a1bSJed Brown int main(int argc, char **argv) 6c4762a1bSJed Brown { 7c4762a1bSJed Brown DM dm; 8c4762a1bSJed Brown PetscSection section; 9c4762a1bSJed Brown PetscFE fe; 1030602db0SMatthew G. Knepley PetscInt dim,c,cStart,cEnd; 11*46181b2aSJed Brown PetscBool view_coord = PETSC_FALSE; 12c4762a1bSJed Brown PetscErrorCode ierr; 13c4762a1bSJed Brown 14c4762a1bSJed Brown ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 15*46181b2aSJed Brown ierr = PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "Tensor closure restrictions", "DMPLEX");CHKERRQ(ierr); 16*46181b2aSJed Brown ierr = PetscOptionsBool("-view_coord", "View coordinates of element closures", "ex8.c", view_coord, &view_coord, NULL);CHKERRQ(ierr); 17*46181b2aSJed Brown ierr = PetscOptionsEnd();CHKERRQ(ierr); 18*46181b2aSJed Brown 1930602db0SMatthew G. Knepley ierr = DMCreate(PETSC_COMM_WORLD, &dm);CHKERRQ(ierr); 2030602db0SMatthew G. Knepley ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 2130602db0SMatthew G. Knepley ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 2230602db0SMatthew G. Knepley ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 2330602db0SMatthew G. Knepley ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 24c4762a1bSJed Brown 25c4762a1bSJed Brown ierr = PetscFECreateDefault(PETSC_COMM_SELF,dim,1,PETSC_FALSE,NULL,PETSC_DETERMINE,&fe);CHKERRQ(ierr); 26c4762a1bSJed Brown ierr = DMAddField(dm,NULL,(PetscObject)fe);CHKERRQ(ierr); 27c4762a1bSJed Brown ierr = DMCreateDS(dm);CHKERRQ(ierr); 28c4762a1bSJed Brown ierr = DMPlexSetClosurePermutationTensor(dm,PETSC_DETERMINE,NULL);CHKERRQ(ierr); 29c4762a1bSJed Brown ierr = DMGetLocalSection(dm,§ion);CHKERRQ(ierr); 30c4762a1bSJed Brown ierr = DMPlexGetHeightStratum(dm,0,&cStart,&cEnd);CHKERRQ(ierr); 31c4762a1bSJed Brown for (c=cStart; c<cEnd; c++) { 32c4762a1bSJed Brown PetscInt numindices,*indices; 3371f0bbf9SMatthew G. Knepley ierr = DMPlexGetClosureIndices(dm,section,section,c,PETSC_TRUE,&numindices,&indices,NULL,NULL);CHKERRQ(ierr); 34*46181b2aSJed Brown ierr = PetscPrintf(PETSC_COMM_SELF,"Element #%" PetscInt_FMT "\n",c-cStart);CHKERRQ(ierr); 35c4762a1bSJed Brown ierr = PetscIntView(numindices,indices,PETSC_VIEWER_STDOUT_SELF);CHKERRQ(ierr); 3671f0bbf9SMatthew G. Knepley ierr = DMPlexRestoreClosureIndices(dm,section,section,c,PETSC_TRUE,&numindices,&indices,NULL,NULL);CHKERRQ(ierr); 37c4762a1bSJed Brown } 38*46181b2aSJed Brown if (view_coord) { 39*46181b2aSJed Brown DM cdm; 40*46181b2aSJed Brown Vec X; 41*46181b2aSJed Brown PetscInt cdim; 42*46181b2aSJed Brown ierr = DMGetCoordinateDM(dm,&cdm);CHKERRQ(ierr); 43*46181b2aSJed Brown ierr = DMPlexSetClosurePermutationTensor(cdm,PETSC_DETERMINE,NULL);CHKERRQ(ierr); 44*46181b2aSJed Brown ierr = DMGetCoordinatesLocal(dm, &X);CHKERRQ(ierr); 45*46181b2aSJed Brown ierr = DMGetCoordinateDim(dm, &cdim);CHKERRQ(ierr); 46*46181b2aSJed Brown for (c=cStart; c<cEnd; c++) { 47*46181b2aSJed Brown PetscScalar *x; 48*46181b2aSJed Brown PetscInt ndof; 49*46181b2aSJed Brown ierr = DMPlexVecGetClosure(cdm, NULL, X, c, &ndof, &x);CHKERRQ(ierr); 50*46181b2aSJed Brown PetscCheck(ndof % cdim == 0, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "ndof not divisible by cdim"); 51*46181b2aSJed Brown ierr = PetscPrintf(PETSC_COMM_SELF,"Element #%" PetscInt_FMT " coordinates\n",c-cStart);CHKERRQ(ierr); 52*46181b2aSJed Brown for (PetscInt i=0; i<ndof; i+= cdim) { 53*46181b2aSJed Brown ierr = PetscScalarView(cdim, &x[i], PETSC_VIEWER_STDOUT_SELF);CHKERRQ(ierr); 54*46181b2aSJed Brown } 55*46181b2aSJed Brown ierr = DMPlexVecRestoreClosure(cdm, NULL, X, c, &ndof, &x);CHKERRQ(ierr); 56*46181b2aSJed Brown } 57*46181b2aSJed Brown } 58c4762a1bSJed Brown ierr = PetscFEDestroy(&fe);CHKERRQ(ierr); 59c4762a1bSJed Brown ierr = DMDestroy(&dm);CHKERRQ(ierr); 60c4762a1bSJed Brown ierr = PetscFinalize(); 61c4762a1bSJed Brown return ierr; 62c4762a1bSJed Brown } 63c4762a1bSJed Brown 64c4762a1bSJed Brown /*TEST 65c4762a1bSJed Brown 66c4762a1bSJed Brown test: 67c4762a1bSJed Brown suffix: 1d_q2 6830602db0SMatthew G. Knepley args: -dm_plex_dim 1 -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_box_faces 2 69c4762a1bSJed Brown test: 70c4762a1bSJed Brown suffix: 2d_q1 7130602db0SMatthew G. Knepley args: -dm_plex_dim 2 -petscspace_degree 1 -dm_plex_simplex 0 -dm_plex_box_faces 2,2 72c4762a1bSJed Brown test: 73c4762a1bSJed Brown suffix: 2d_q2 7430602db0SMatthew G. Knepley args: -dm_plex_dim 2 -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_box_faces 2,2 75c4762a1bSJed Brown test: 76c4762a1bSJed Brown suffix: 2d_q3 7730602db0SMatthew G. Knepley args: -dm_plex_dim 2 -petscspace_degree 3 -dm_plex_simplex 0 -dm_plex_box_faces 1,1 78c4762a1bSJed Brown test: 79c4762a1bSJed Brown suffix: 3d_q1 8030602db0SMatthew G. Knepley args: -dm_plex_dim 3 -petscspace_degree 1 -dm_plex_simplex 0 -dm_plex_box_faces 1,1,1 81*46181b2aSJed Brown test: 82*46181b2aSJed Brown suffix: 1d_q1_periodic 83*46181b2aSJed Brown args: -dm_plex_dim 1 -petscspace_degree 1 -dm_plex_simplex 0 -dm_plex_box_faces 2 -dm_plex_box_bd periodic -dm_view 84*46181b2aSJed Brown test: 85*46181b2aSJed Brown suffix: 2d_q1_periodic 86*46181b2aSJed Brown args: -dm_plex_dim 2 -petscspace_degree 1 -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -dm_plex_box_bd periodic,none -dm_view 87*46181b2aSJed Brown test: 88*46181b2aSJed Brown suffix: 3d_q2_periodic 89*46181b2aSJed Brown args: -dm_plex_dim 3 -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_box_faces 2,2,2 -dm_plex_box_bd periodic,none,periodic -dm_view 90c4762a1bSJed Brown 91c4762a1bSJed Brown TEST*/ 92