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