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