xref: /petsc/src/dm/impls/plex/tests/ex13.c (revision f97672e55eacc8688507b9471cd7ec2664d7f203)
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 
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(0);
22 }
23 
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         case 1: {
57           PetscInt triSizes_p2[2]  = {6, 2};
58           PetscInt triPoints_p2[8] = {1, 2, 3, 4, 6, 7, 0, 5};
59 
60           PetscCall(PetscMalloc2(2, &sizes, 8, &points));
61           PetscCall(PetscArraycpy(sizes,  triSizes_p2, 2));
62           PetscCall(PetscArraycpy(points, triPoints_p2, 8));
63           break;}
64         default:
65           SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for triangular mesh on 2 procs", user->testNum);
66         }
67       } else if (dim == 2 && simplex && size == 3) {
68         PetscInt triSizes_p3[3]  = {3, 3, 2};
69         PetscInt triPoints_p3[8] = {1, 2, 4, 3, 6, 7, 0, 5};
70 
71         PetscCall(PetscMalloc2(3, &sizes, 8, &points));
72         PetscCall(PetscArraycpy(sizes,  triSizes_p3, 3));
73         PetscCall(PetscArraycpy(points, triPoints_p3, 8));
74       } else if (dim == 2 && !simplex && size == 2) {
75         PetscInt quadSizes_p2[2]  = {2, 2};
76         PetscInt quadPoints_p2[4] = {2, 3, 0, 1};
77 
78         PetscCall(PetscMalloc2(2, &sizes, 4, &points));
79         PetscCall(PetscArraycpy(sizes,  quadSizes_p2, 2));
80         PetscCall(PetscArraycpy(points, quadPoints_p2, 4));
81       } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition");
82     }
83     PetscCall(DMPlexGetPartitioner(*dm, &part));
84     PetscCall(PetscPartitionerSetType(part, PETSCPARTITIONERSHELL));
85     PetscCall(PetscPartitionerShellSetPartition(part, size, sizes, points));
86     PetscCall(PetscFree2(sizes, points));
87   }
88   PetscCall(DMPlexDistribute(*dm, 0, NULL, &dmDist));
89   if (dmDist) {
90     PetscCall(DMDestroy(dm));
91     *dm  = dmDist;
92   }
93   PetscCall(PetscObjectSetName((PetscObject) *dm, simplex ? "Simplicial Mesh" : "Tensor Product Mesh"));
94   PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
95   PetscFunctionReturn(0);
96 }
97 
98 static PetscErrorCode ScrambleOrientation(DM dm, AppCtx *user)
99 {
100   PetscInt       h, cStart, cEnd, c;
101 
102   PetscFunctionBeginUser;
103   PetscCall(DMPlexGetVTKCellHeight(dm, &h));
104   PetscCall(DMPlexGetHeightStratum(dm, h, &cStart, &cEnd));
105   for (c = cStart; c < cEnd; ++c) {
106     /* Could use PetscRand instead */
107     if (c%2) PetscCall(DMPlexOrientPoint(dm, c, -1));
108   }
109   PetscFunctionReturn(0);
110 }
111 
112 static PetscErrorCode TestOrientation(DM dm, AppCtx *user)
113 {
114   PetscFunctionBeginUser;
115   PetscCall(ScrambleOrientation(dm, user));
116   PetscCall(DMPlexOrient(dm));
117   PetscCall(DMViewFromOptions(dm, NULL, "-oriented_dm_view"));
118   PetscFunctionReturn(0);
119 }
120 
121 int main(int argc, char **argv)
122 {
123   DM             dm;
124   AppCtx         user; /* user-defined work context */
125 
126   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
127   PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
128   PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
129   PetscCall(TestOrientation(dm, &user));
130   PetscCall(DMDestroy(&dm));
131   PetscCall(PetscFinalize());
132   return 0;
133 }
134 
135 /*TEST
136   testset:
137     requires: triangle
138     args: -dm_coord_space 0 -dm_view ascii::ascii_info_detail -oriented_dm_view ascii::ascii_info_detail -orientation_view
139 
140     test:
141       suffix: 0
142       args: -test_partition 0
143     test:
144       suffix: 1
145       nsize: 2
146     test:
147       suffix: 2
148       nsize: 2
149       args: -test_num 1
150     test:
151       suffix: 3
152       nsize: 3
153       args: -orientation_view_synchronized
154 
155 TEST*/
156