Actual source code: ex8.c

  1: static char help[] = "Element closure restrictions in tensor/lexicographic/spectral-element ordering using DMPlex\n\n";

  3: #include <petscdmplex.h>

  5: int main(int argc, char **argv)
  6: {
  7:   DM             dm;
  8:   PetscSection   section;
  9:   PetscFE        fe;
 10:   PetscInt       dim,c,cStart,cEnd;

 13:   PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr;
 14:   DMCreate(PETSC_COMM_WORLD, &dm);
 15:   DMSetType(dm, DMPLEX);
 16:   DMSetFromOptions(dm);
 17:   DMViewFromOptions(dm, NULL, "-dm_view");
 18:   DMGetDimension(dm, &dim);

 20:   PetscFECreateDefault(PETSC_COMM_SELF,dim,1,PETSC_FALSE,NULL,PETSC_DETERMINE,&fe);
 21:   DMSetFromOptions(dm);
 22:   DMAddField(dm,NULL,(PetscObject)fe);
 23:   DMCreateDS(dm);
 24:   DMPlexSetClosurePermutationTensor(dm,PETSC_DETERMINE,NULL);
 25:   DMGetLocalSection(dm,&section);
 26:   DMPlexGetHeightStratum(dm,0,&cStart,&cEnd);
 27:   for (c=cStart; c<cEnd; c++) {
 28:     PetscInt numindices,*indices;
 29:     DMPlexGetClosureIndices(dm,section,section,c,PETSC_TRUE,&numindices,&indices,NULL,NULL);
 30:     PetscPrintf(PETSC_COMM_SELF,"Element #%D\n",c-cStart);
 31:     PetscIntView(numindices,indices,PETSC_VIEWER_STDOUT_SELF);
 32:     DMPlexRestoreClosureIndices(dm,section,section,c,PETSC_TRUE,&numindices,&indices,NULL,NULL);
 33:   }

 35:   PetscFEDestroy(&fe);
 36:   DMDestroy(&dm);
 37:   PetscFinalize();
 38:   return ierr;
 39: }

 41: /*TEST

 43:   test:
 44:     suffix: 1d_q2
 45:     args: -dm_plex_dim 1 -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_box_faces 2
 46:   test:
 47:     suffix: 2d_q1
 48:     args: -dm_plex_dim 2 -petscspace_degree 1 -dm_plex_simplex 0 -dm_plex_box_faces 2,2
 49:   test:
 50:     suffix: 2d_q2
 51:     args: -dm_plex_dim 2 -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_box_faces 2,2
 52:   test:
 53:     suffix: 2d_q3
 54:     args: -dm_plex_dim 2 -petscspace_degree 3 -dm_plex_simplex 0 -dm_plex_box_faces 1,1
 55:   test:
 56:     suffix: 3d_q1
 57:     args: -dm_plex_dim 3 -petscspace_degree 1 -dm_plex_simplex 0 -dm_plex_box_faces 1,1,1

 59: TEST*/