Actual source code: ex13.c
1: static char help[] = "Orient a mesh in parallel\n\n";
3: #include <petscdmplex.h>
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;
11: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
12: {
13: PetscFunctionBeginUser;
14: options->testPartition = PETSC_TRUE;
15: options->testNum = 0;
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: }
24: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
25: {
26: DM dmDist = NULL;
27: PetscBool simplex;
28: PetscInt dim;
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;
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};
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};
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};
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};
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: }
100: static PetscErrorCode ScrambleOrientation(DM dm, AppCtx *user)
101: {
102: PetscInt h, cStart, cEnd, c;
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: }
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: }
123: int main(int argc, char **argv)
124: {
125: DM dm;
126: AppCtx user; /* user-defined work context */
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: }
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
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
158: TEST*/