1 static char help[] = "Tests for submesh creation for periodic meshes\n\n"; 2 3 #include <petscdmplex.h> 4 #include <petsc/private/dmimpl.h> 5 #include <petscsf.h> 6 7 /* 3 x 2 2D mesh periodic in x-direction using nproc = 1. 8 9 Cell submesh (subdim = 2): 10 11 12--21--13--22--14--23---~ 6--12---7--13---8--14---~ 12 | | | ~ | | | ~ 13 25 3 27 4 29 5 ~ 15 0 16 1 17 2 ~ 14 | | | ~ | | | ~ 15 9--18--10--19--11--20---~ -> 3---9---4--10---5--11---~ 16 | | | ~ 17 24 0 26 1 28 2 ~ 18 | | | ~ 19 6--15---7--16---8--17---~ 20 21 Facet submesh (subdim = 1): 22 23 12--21--13--22--14--23---~ 3---0---4---1---5---2---~ 24 | | | ~ 25 25 3 27 4 29 5 ~ 26 | | | ~ 27 9--18--10--19--11--20---~ -> 28 | | | ~ 29 24 0 26 1 28 2 ~ 30 | | | ~ 31 6--15---7--16---8--17---~ 32 33 */ 34 35 typedef struct { 36 PetscInt subdim; /* The topological dimension of the submesh */ 37 } AppCtx; 38 39 PetscErrorCode ProcessOptions(AppCtx *options) 40 { 41 PetscFunctionBegin; 42 options->subdim = 2; 43 44 PetscOptionsBegin(PETSC_COMM_SELF, "", "Filtering Problem Options", "DMPLEX"); 45 PetscCall(PetscOptionsBoundedInt("-subdim", "The topological dimension of the submesh", "ex74.c", options->subdim, &options->subdim, NULL, 0)); 46 PetscOptionsEnd(); 47 PetscFunctionReturn(PETSC_SUCCESS); 48 } 49 50 int main(int argc, char **argv) 51 { 52 PetscInt dim = 2; 53 DM dm, subdm, coorddm; 54 PetscFE coordfe; 55 const PetscInt faces[2] = {3, 2}; 56 const PetscReal lower[2] = {0., 0.}, upper[2] = {3., 2.}; 57 const DMBoundaryType periodicity[2] = {DM_BOUNDARY_PERIODIC, DM_BOUNDARY_NONE}; 58 DMLabel filter; 59 const PetscInt filterValue = 1; 60 MPI_Comm comm; 61 PetscMPIInt size; 62 AppCtx user; 63 64 PetscFunctionBeginUser; 65 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 66 PetscCall(ProcessOptions(&user)); 67 comm = PETSC_COMM_WORLD; 68 PetscCallMPI(MPI_Comm_size(comm, &size)); 69 if (size != 1) { 70 PetscCall(PetscPrintf(comm, "This example is specifically designed for comm size == 1.\n")); 71 PetscCall(PetscFinalize()); 72 return 0; 73 } 74 PetscCall(DMPlexCreateBoxMesh(comm, dim, PETSC_FALSE, faces, lower, upper, periodicity, PETSC_TRUE, 0, PETSC_TRUE, &dm)); 75 /* Reset DG coordinates so that DMLocalizeCoordinates() will run again */ 76 /* in DMSetFromOptions() after CG coordinate FE is set. */ 77 dm->coordinates[1].dim = -1; 78 /* Localize coordinates on facets, too, if we are to make a facet submesh. */ 79 /* Otherwise, DG coordinate values will not be copied from the parent mesh. */ 80 PetscCall(DMSetFromOptions(dm)); 81 PetscCall(DMLabelCreate(PETSC_COMM_SELF, "filter", &filter)); 82 switch (user.subdim) { 83 case 2: 84 PetscCall(DMLabelSetValue(filter, 3, filterValue)); 85 PetscCall(DMLabelSetValue(filter, 4, filterValue)); 86 PetscCall(DMLabelSetValue(filter, 5, filterValue)); 87 break; 88 case 1: 89 PetscCall(DMLabelSetValue(filter, 21, filterValue)); 90 PetscCall(DMLabelSetValue(filter, 22, filterValue)); 91 PetscCall(DMLabelSetValue(filter, 23, filterValue)); 92 break; 93 default: 94 PetscCall(PetscPrintf(comm, "This example is only for subdim == {1, 2}.\n")); 95 PetscCall(PetscFinalize()); 96 return 0; 97 } 98 PetscCall(DMPlexFilter(dm, filter, filterValue, PETSC_FALSE, PETSC_FALSE, NULL, &subdm)); 99 PetscCall(DMLabelDestroy(&filter)); 100 PetscCall(DMDestroy(&dm)); 101 PetscCall(PetscObjectSetName((PetscObject)subdm, "Example_SubDM")); 102 PetscCall(DMViewFromOptions(subdm, NULL, "-dm_view")); 103 PetscCall(DMGetCellCoordinateDM(subdm, &coorddm)); 104 PetscCall(DMGetField(coorddm, 0, NULL, (PetscObject *)&coordfe)); 105 PetscCall(PetscFEView(coordfe, PETSC_VIEWER_STDOUT_WORLD)); 106 PetscCall(DMDestroy(&subdm)); 107 PetscCall(PetscFinalize()); 108 return 0; 109 } 110 111 /*TEST 112 113 testset: 114 args: -dm_coord_space 1 -dm_coord_petscspace_degree 1 -dm_localize 1 -dm_view ascii::ascii_info_detail 115 116 test: 117 suffix: 0 118 args: -dm_sparse_localize 0 -dm_localize_height 0 -subdim 2 119 120 test: 121 suffix: 1 122 args: -dm_sparse_localize 1 -dm_localize_height 0 -subdim 2 123 124 test: 125 suffix: 2 126 args: -dm_sparse_localize 0 -dm_localize_height 1 -subdim 1 127 128 test: 129 suffix: 3 130 args: -dm_sparse_localize 1 -dm_localize_height 1 -subdim 1 131 132 TEST*/ 133