1 static char help[] = "Orient a mesh in parallel\n\n";
2
3 #include <petscdmplex.h>
4
5 typedef struct {
6 /* Domain and mesh definition */
7 PetscInt dim; /* The topological mesh dimension */
8 PetscBool cellSimplex; /* Use simplices or hexes */
9 char filename[PETSC_MAX_PATH_LEN]; /* Import mesh from file */
10 PetscBool testPartition; /* Use a fixed partitioning for testing */
11 PetscInt testNum; /* Labels the different test partitions */
12 } AppCtx;
13
ProcessOptions(MPI_Comm comm,AppCtx * options)14 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
15 {
16 PetscErrorCode ierr;
17
18 PetscFunctionBeginUser;
19 options->dim = 2;
20 options->cellSimplex = PETSC_TRUE;
21 options->filename[0] = '\0';
22 options->testPartition = PETSC_TRUE;
23 options->testNum = 0;
24
25 ierr = PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");CHKERRQ(ierr);
26 ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex13.c", options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr);
27 ierr = PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex13.c", options->cellSimplex, &options->cellSimplex, NULL);CHKERRQ(ierr);
28 ierr = PetscOptionsString("-filename", "The mesh file", "ex13.c", options->filename, options->filename, sizeof(options->filename), NULL);CHKERRQ(ierr);
29 ierr = PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex13.c", options->testPartition, &options->testPartition, NULL);CHKERRQ(ierr);
30 ierr = PetscOptionsBoundedInt("-test_num", "The test partition number", "ex13.c", options->testNum, &options->testNum, NULL,0);CHKERRQ(ierr);
31 ierr = PetscOptionsEnd();
32 PetscFunctionReturn(0);
33 }
34
CreateMesh(MPI_Comm comm,AppCtx * user,DM * dm)35 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
36 {
37 DM dmDist = NULL;
38 PetscInt dim = user->dim;
39 PetscBool cellSimplex = user->cellSimplex;
40 const char *filename = user->filename;
41 size_t len;
42 PetscErrorCode ierr;
43
44 PetscFunctionBeginUser;
45 ierr = PetscStrlen(filename, &len);CHKERRQ(ierr);
46 if (len) {ierr = DMPlexCreateFromFile(comm, filename, PETSC_TRUE, dm);CHKERRQ(ierr);}
47 else {ierr = DMPlexCreateBoxMesh(comm, dim, cellSimplex, NULL, NULL, NULL, NULL, PETSC_TRUE, dm);CHKERRQ(ierr);}
48 if (user->testPartition) {
49 PetscPartitioner part;
50 PetscInt *sizes = NULL;
51 PetscInt *points = NULL;
52 PetscMPIInt rank, size;
53
54 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr);
55 ierr = MPI_Comm_size(comm, &size);CHKERRQ(ierr);
56 if (!rank) {
57 if (dim == 2 && cellSimplex && size == 2) {
58 switch (user->testNum) {
59 case 0: {
60 PetscInt triSizes_p2[2] = {4, 4};
61 PetscInt triPoints_p2[8] = {3, 5, 6, 7, 0, 1, 2, 4};
62
63 ierr = PetscMalloc2(2, &sizes, 8, &points);CHKERRQ(ierr);
64 ierr = PetscArraycpy(sizes, triSizes_p2, 2);CHKERRQ(ierr);
65 ierr = PetscArraycpy(points, triPoints_p2, 8);CHKERRQ(ierr);
66 break;}
67 case 1: {
68 PetscInt triSizes_p2[2] = {6, 2};
69 PetscInt triPoints_p2[8] = {1, 2, 3, 4, 6, 7, 0, 5};
70
71 ierr = PetscMalloc2(2, &sizes, 8, &points);CHKERRQ(ierr);
72 ierr = PetscArraycpy(sizes, triSizes_p2, 2);CHKERRQ(ierr);
73 ierr = PetscArraycpy(points, triPoints_p2, 8);CHKERRQ(ierr);
74 break;}
75 default:
76 SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for triangular mesh on 2 procs", user->testNum);
77 }
78 } else if (dim == 2 && cellSimplex && size == 3) {
79 PetscInt triSizes_p3[3] = {3, 3, 2};
80 PetscInt triPoints_p3[8] = {1, 2, 4, 3, 6, 7, 0, 5};
81
82 ierr = PetscMalloc2(3, &sizes, 8, &points);CHKERRQ(ierr);
83 ierr = PetscArraycpy(sizes, triSizes_p3, 3);CHKERRQ(ierr);
84 ierr = PetscArraycpy(points, triPoints_p3, 8);CHKERRQ(ierr);
85 } else if (dim == 2 && !cellSimplex && size == 2) {
86 PetscInt quadSizes_p2[2] = {2, 2};
87 PetscInt quadPoints_p2[4] = {2, 3, 0, 1};
88
89 ierr = PetscMalloc2(2, &sizes, 4, &points);CHKERRQ(ierr);
90 ierr = PetscArraycpy(sizes, quadSizes_p2, 2);CHKERRQ(ierr);
91 ierr = PetscArraycpy(points, quadPoints_p2, 4);CHKERRQ(ierr);
92 } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition");
93 }
94 ierr = DMPlexGetPartitioner(*dm, &part);CHKERRQ(ierr);
95 ierr = PetscPartitionerSetType(part, PETSCPARTITIONERSHELL);CHKERRQ(ierr);
96 ierr = PetscPartitionerShellSetPartition(part, size, sizes, points);CHKERRQ(ierr);
97 ierr = PetscFree2(sizes, points);CHKERRQ(ierr);
98 }
99 ierr = DMPlexDistribute(*dm, 0, NULL, &dmDist);CHKERRQ(ierr);
100 if (dmDist) {
101 ierr = DMDestroy(dm);CHKERRQ(ierr);
102 *dm = dmDist;
103 }
104 ierr = PetscObjectSetName((PetscObject) *dm, cellSimplex ? "Simplicial Mesh" : "Tensor Product Mesh");CHKERRQ(ierr);
105 ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr);
106 PetscFunctionReturn(0);
107 }
108
ScrambleOrientation(DM dm,AppCtx * user)109 static PetscErrorCode ScrambleOrientation(DM dm, AppCtx *user)
110 {
111 PetscInt h, cStart, cEnd, c;
112 PetscErrorCode ierr;
113
114 PetscFunctionBeginUser;
115 ierr = DMPlexGetVTKCellHeight(dm, &h);CHKERRQ(ierr);
116 ierr = DMPlexGetHeightStratum(dm, h, &cStart, &cEnd);CHKERRQ(ierr);
117 for (c = cStart; c < cEnd; ++c) {
118 /* Could use PetscRand instead */
119 if (c%2) {ierr = DMPlexReverseCell(dm, c);CHKERRQ(ierr);}
120 }
121 PetscFunctionReturn(0);
122 }
123
TestOrientation(DM dm,AppCtx * user)124 static PetscErrorCode TestOrientation(DM dm, AppCtx *user)
125 {
126 PetscErrorCode ierr;
127
128 PetscFunctionBeginUser;
129 ierr = ScrambleOrientation(dm, user);CHKERRQ(ierr);
130 ierr = DMPlexOrient(dm);CHKERRQ(ierr);
131 ierr = DMViewFromOptions(dm, NULL, "-oriented_dm_view");CHKERRQ(ierr);
132 PetscFunctionReturn(0);
133 }
134
main(int argc,char ** argv)135 int main(int argc, char **argv)
136 {
137 DM dm;
138 AppCtx user; /* user-defined work context */
139 PetscErrorCode ierr;
140
141 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr;
142 ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr);
143 ierr = CreateMesh(PETSC_COMM_WORLD, &user, &dm);CHKERRQ(ierr);
144 ierr = TestOrientation(dm, &user);CHKERRQ(ierr);
145 ierr = DMDestroy(&dm);CHKERRQ(ierr);
146 ierr = PetscFinalize();
147 return ierr;
148 }
149
150 /*TEST
151 test:
152 suffix: 0
153 requires: triangle
154 args: -test_partition 0 -dm_view ascii::ascii_info_detail -oriented_dm_view ascii::ascii_info_detail -orientation_view
155 test:
156 suffix: 1
157 requires: triangle
158 nsize: 2
159 args: -dm_view ascii::ascii_info_detail -oriented_dm_view ascii::ascii_info_detail -orientation_view
160 test:
161 suffix: 2
162 requires: triangle
163 nsize: 2
164 args: -test_num 1 -dm_view ascii::ascii_info_detail -oriented_dm_view ascii::ascii_info_detail -orientation_view
165 test:
166 suffix: 3
167 requires: triangle
168 nsize: 3
169 args: -dm_view ascii::ascii_info_detail -oriented_dm_view ascii::ascii_info_detail -orientation_view -orientation_view_synchronized
170
171 TEST*/
172