Actual source code: ex71.c
1: static char help[] = "Tests for submesh creation\n\n";
3: #include <petscdmplex.h>
4: #include <petscsf.h>
6: /* Submesh of a 2 x 2 mesh using 3 processes.
8: Local numbering on each rank:
10: (6)(16)-(7)(17)-(8) 5--14---6--15---7 (10)(20)(11)(21)(12)
11: | | | | | | | | |
12: (18) (1)(19) (2)(20) 16 0 17 1 18 (22) (2)(23) (3)(24)
13: | | | | | | | | |
14: (5)(15)(11)(22)(12) 4--13-(11)(22)(12) (9)(19)--6--14---7
15: | | | | | | | | |
16: 14 0 (23) (3)(24) (20) (2)(23) (3)(24) (18) (1) 15 0 16
17: | | | | | | | | |
18: 4--13--(9)(21)(10) (8)(19)-(9)(21)(10) (8)(17)--4--13---5
20: mesh_0 mesh_1 mesh_2
22: where () represents ghost points. We extract the left 2 cells.
23: With sanitize_submesh = PETSC_FALSE, we get:
25: (4)(11)-(5) 3---9---4
26: | | | |
27: (12) (1)(13) 10 0 11
28: | | | |
29: (3)(10) (7) 2---8---7
30: | | | |
31: 9 0 (14) (13) (1) 14
32: | | | |
33: 2---8--(6) (5)(12)--6
35: On the other hand, with sanitize_submesh = PETSC_TRUE, we get:
37: (4)(11)-(5) 3---9---4
38: | | | |
39: (12) (1)(13) 10 0 11
40: | | | |
41: (3)(10) (7) 2---8---7
42: | | | |
43: 9 0 14 (13) (1)(14)
44: | | | |
45: 2---8---6 (5)(12)-(6)
47: submesh_0 submesh_1 submesh_2
49: as points 15 and 4 of mesh_2 are in the closure of a submesh cell owned by rank 0 (point 0 of submesh_0),
50: and not in the closure of any submesh cell owned by rank 1.
52: */
54: typedef struct {
55: PetscBool ignoreLabelHalo; /* Ignore filter values in the halo. */
56: PetscBool sanitizeSubmesh; /* Sanitize submesh. */
57: } AppCtx;
59: PetscErrorCode ProcessOptions(AppCtx *options)
60: {
61: PetscFunctionBegin;
62: options->ignoreLabelHalo = PETSC_FALSE;
63: options->sanitizeSubmesh = PETSC_FALSE;
65: PetscOptionsBegin(PETSC_COMM_SELF, "", "Filtering Problem Options", "DMPLEX");
66: PetscCall(PetscOptionsBool("-ignore_label_halo", "Ignore filter values in the halo", "ex80.c", options->ignoreLabelHalo, &options->ignoreLabelHalo, NULL));
67: PetscCall(PetscOptionsBool("-sanitize_submesh", "Sanitize submesh", "ex80.c", options->sanitizeSubmesh, &options->sanitizeSubmesh, NULL));
68: PetscOptionsEnd();
69: PetscFunctionReturn(PETSC_SUCCESS);
70: }
71: int main(int argc, char **argv)
72: {
73: DM dm, subdm;
74: PetscSF ownershipTransferSF;
75: DMLabel filter;
76: const PetscInt filterValue = 1;
77: MPI_Comm comm;
78: PetscMPIInt size, rank;
79: AppCtx user;
81: PetscFunctionBeginUser;
82: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
83: PetscCall(ProcessOptions(&user));
84: comm = PETSC_COMM_WORLD;
85: PetscCallMPI(MPI_Comm_size(comm, &size));
86: if (size != 3) {
87: PetscCall(PetscPrintf(comm, "This example is specifically designed for comm size == 3.\n"));
88: PetscCall(PetscFinalize());
89: return 0;
90: }
91: PetscCallMPI(MPI_Comm_rank(comm, &rank));
92: {
93: DM pdm;
94: const PetscInt faces[2] = {2, 2};
95: PetscInt overlap = 1;
97: PetscCall(DMPlexCreateBoxMesh(comm, 2, PETSC_FALSE, faces, NULL, NULL, NULL, PETSC_TRUE, 0, PETSC_TRUE, &dm));
98: {
99: PetscPartitioner part;
100: PetscInt *sizes = NULL;
101: PetscInt *points = NULL;
103: if (rank == 0) {
104: PetscInt sizes1[3] = {1, 2, 1};
105: PetscInt points1[4] = {0, 2, 3, 1};
107: PetscCall(PetscMalloc2(3, &sizes, 4, &points));
108: PetscCall(PetscArraycpy(sizes, sizes1, 3));
109: PetscCall(PetscArraycpy(points, points1, 4));
110: }
111: PetscCall(DMPlexGetPartitioner(dm, &part));
112: PetscCall(PetscPartitionerSetType(part, PETSCPARTITIONERSHELL));
113: PetscCall(PetscPartitionerShellSetPartition(part, size, sizes, points));
114: PetscCall(PetscFree2(sizes, points));
115: }
116: PetscCall(DMSetAdjacency(dm, -1, PETSC_FALSE, PETSC_TRUE));
117: PetscCall(DMPlexDistribute(dm, overlap, NULL, &pdm));
118: if (pdm) {
119: PetscCall(DMDestroy(&dm));
120: dm = pdm;
121: }
122: }
123: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "filter", &filter));
124: switch (rank) {
125: case 0:
126: PetscCall(DMLabelSetValue(filter, 0, filterValue));
127: PetscCall(DMLabelSetValue(filter, 1, filterValue));
128: break;
129: case 1:
130: PetscCall(DMLabelSetValue(filter, 0, filterValue));
131: PetscCall(DMLabelSetValue(filter, 2, filterValue));
132: break;
133: case 2:
134: break;
135: }
136: PetscCall(PetscObjectSetName((PetscObject)dm, "Example_DM"));
137: PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
138: PetscCall(DMPlexFilter(dm, filter, filterValue, user.ignoreLabelHalo, user.sanitizeSubmesh, &ownershipTransferSF, &subdm));
139: PetscCall(DMLabelDestroy(&filter));
140: PetscCall(DMDestroy(&dm));
141: PetscCall(PetscObjectSetName((PetscObject)subdm, "Example_SubDM"));
142: PetscCall(DMViewFromOptions(subdm, NULL, "-dm_view"));
143: PetscCall(DMDestroy(&subdm));
144: PetscCall(PetscObjectSetName((PetscObject)ownershipTransferSF, "Example_Ownership_Transfer_SF"));
145: PetscCall(PetscSFView(ownershipTransferSF, PETSC_VIEWER_STDOUT_WORLD));
146: PetscCall(PetscSFDestroy(&ownershipTransferSF));
147: PetscCall(PetscFinalize());
148: return 0;
149: }
151: /*TEST
153: testset:
154: nsize: 3
155: args: -dm_view ascii::ascii_info_detail
157: test:
158: suffix: 0
159: args:
161: test:
162: suffix: 1
163: args: -sanitize_submesh
165: test:
166: suffix: 2
167: args: -ignore_label_halo
169: TEST*/