Actual source code: ex44.c
1: static const char help[] = "Tests for mesh extrusion";
3: #include <petscdmplex.h>
4: #include <petscdmforest.h>
6: typedef struct {
7: char bdLabel[PETSC_MAX_PATH_LEN]; /* The boundary label name */
8: PetscInt Nbd; /* The number of boundary markers to extrude, 0 for all */
9: PetscInt bd[64]; /* The boundary markers to be extruded */
10: PetscBool testForest; /* Convert the mesh to a forest and back */
11: } AppCtx;
13: PETSC_EXTERN PetscErrorCode pyramidNormal(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *);
15: /* The pyramid apex is at (0.5, 0.5, -1) */
16: PetscErrorCode pyramidNormal(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt r, PetscScalar u[], PetscCtx ctx)
17: {
18: PetscReal apex[3] = {0.5, 0.5, -1.0};
19: PetscInt d;
21: for (d = 0; d < dim; ++d) u[d] = x[d] - apex[d];
22: for (d = dim; d < 3; ++d) u[d] = 0.0 - apex[d];
23: return PETSC_SUCCESS;
24: }
26: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
27: {
28: PetscInt n = 64;
29: PetscBool flg;
31: PetscFunctionBeginUser;
32: PetscCall(PetscStrncpy(options->bdLabel, "marker", sizeof(options->bdLabel)));
33: PetscOptionsBegin(comm, "", "Parallel Mesh Adaptation Options", "DMPLEX");
34: PetscCall(PetscOptionsString("-label", "The boundary label name", "ex44.c", options->bdLabel, options->bdLabel, sizeof(options->bdLabel), NULL));
35: PetscCall(PetscOptionsIntArray("-bd", "The boundaries to be extruded", "ex44.c", options->bd, &n, &flg));
36: options->Nbd = flg ? n : 0;
37: options->testForest = PETSC_FALSE;
38: PetscCall(PetscOptionsBool("-test_forest", "Create a DMForest and then convert to DMPlex before testing", "ex44.c", PETSC_FALSE, &options->testForest, NULL));
39: PetscOptionsEnd();
40: PetscFunctionReturn(PETSC_SUCCESS);
41: }
43: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *ctx, DM *dm)
44: {
45: PetscFunctionBegin;
46: PetscCall(DMCreate(comm, dm));
47: PetscCall(DMSetType(*dm, DMPLEX));
48: PetscCall(DMSetFromOptions(*dm));
49: if (ctx->testForest) {
50: DM dmForest;
51: PetscInt dim;
52: PetscCall(DMGetDimension(*dm, &dim));
53: PetscCall(DMCreate(PETSC_COMM_WORLD, &dmForest));
54: PetscCall(DMSetType(dmForest, (dim == 2) ? DMP4EST : DMP8EST));
55: PetscCall(DMForestSetBaseDM(dmForest, *dm));
56: PetscCall(DMSetUp(dmForest));
57: PetscCall(DMDestroy(dm));
58: PetscCall(DMConvert(dmForest, DMPLEX, dm));
59: PetscCall(DMDestroy(&dmForest));
60: PetscCall(PetscObjectSetName((PetscObject)*dm, "forest"));
61: }
62: PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
63: PetscFunctionReturn(PETSC_SUCCESS);
64: }
66: static PetscErrorCode CreateAdaptLabel(DM dm, AppCtx *ctx, DMLabel *adaptLabel)
67: {
68: DMLabel label;
69: PetscInt b;
71: PetscFunctionBegin;
72: if (!ctx->Nbd) {
73: *adaptLabel = NULL;
74: PetscFunctionReturn(PETSC_SUCCESS);
75: }
76: PetscCall(DMGetLabel(dm, ctx->bdLabel, &label));
77: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Adaptation Label", adaptLabel));
78: for (b = 0; b < ctx->Nbd; ++b) {
79: IS bdIS;
80: const PetscInt *points;
81: PetscInt n, i;
83: PetscCall(DMLabelGetStratumIS(label, ctx->bd[b], &bdIS));
84: if (!bdIS) continue;
85: PetscCall(ISGetLocalSize(bdIS, &n));
86: PetscCall(ISGetIndices(bdIS, &points));
87: for (i = 0; i < n; ++i) PetscCall(DMLabelSetValue(*adaptLabel, points[i], DM_ADAPT_REFINE));
88: PetscCall(ISRestoreIndices(bdIS, &points));
89: PetscCall(ISDestroy(&bdIS));
90: }
91: PetscFunctionReturn(PETSC_SUCCESS);
92: }
94: int main(int argc, char **argv)
95: {
96: DM dm, dma;
97: DMLabel adaptLabel;
98: AppCtx ctx;
100: PetscFunctionBeginUser;
101: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
102: PetscCall(ProcessOptions(PETSC_COMM_WORLD, &ctx));
103: PetscCall(CreateMesh(PETSC_COMM_WORLD, &ctx, &dm));
104: PetscCall(CreateAdaptLabel(dm, &ctx, &adaptLabel));
105: if (adaptLabel) PetscCall(DMAdaptLabel(dm, adaptLabel, &dma));
106: else PetscCall(DMExtrude(dm, 3, &dma));
107: PetscCall(PetscObjectSetName((PetscObject)dma, "Adapted Mesh"));
108: PetscCall(DMLabelDestroy(&adaptLabel));
109: PetscCall(DMDestroy(&dm));
110: PetscCall(DMViewFromOptions(dma, NULL, "-adapt_dm_view"));
111: PetscCall(DMDestroy(&dma));
112: PetscCall(PetscFinalize());
113: return 0;
114: }
116: /*TEST
118: testset:
119: args: -dm_view -adapt_dm_view -dm_plex_check_all
121: test:
122: suffix: seg_periodic_0
123: args: -dm_plex_dim 1 -dm_plex_box_faces 3 -dm_plex_transform_extrude_periodic -dm_plex_transform_extrude_use_tensor 0
125: test:
126: suffix: seg_periodic_1
127: args: -dm_plex_dim 1 -dm_plex_box_faces 4 -dm_plex_box_bd periodic \
128: -dm_plex_transform_extrude_layers 3 -dm_plex_transform_extrude_use_tensor 0
130: test:
131: suffix: tri_tensor_0
132: requires: triangle
133: args: -dm_plex_transform_extrude_use_tensor {{0 1}separate output}
135: test:
136: suffix: quad_tensor_0
137: args: -dm_plex_simplex 0 -dm_plex_transform_extrude_use_tensor {{0 1}separate output}
139: test:
140: suffix: quad_tensor_0_forest
141: requires: p4est
142: args: -dm_plex_simplex 0 -dm_plex_transform_extrude_use_tensor {{0 1}separate output} -test_forest 1
144: test:
145: suffix: quad_normal_0
146: args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal 0,1,1
148: test:
149: suffix: quad_normal_0_forest
150: requires: p4est
151: args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal 0,1,1 -test_forest 1
153: test:
154: suffix: quad_normal_1
155: args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal_function pyramidNormal
157: test:
158: suffix: quad_normal_1_forest
159: requires: p4est
160: args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal_function pyramidNormal -test_forest 1
162: test:
163: suffix: quad_symmetric_0
164: args: -dm_plex_simplex 0 -dm_plex_transform_extrude_symmetric
166: test:
167: suffix: quad_symmetric_0_forest
168: requires: p4est
169: args: -dm_plex_simplex 0 -dm_plex_transform_extrude_symmetric -test_forest 1
171: test:
172: suffix: quad_label
173: args: -dm_plex_simplex 0 -dm_plex_transform_label_replica_inc {{0 100}separate output}
175: test:
176: suffix: quad_periodic_0
177: args: -dm_plex_simplex 0 -dm_plex_transform_extrude_periodic -dm_plex_transform_extrude_use_tensor 0
179: testset:
180: args: -dm_adaptor cellrefiner -dm_plex_transform_type extrude \
181: -dm_view -adapt_dm_view -dm_plex_check_all
183: test:
184: suffix: tri_adapt_0
185: requires: triangle
186: args: -dm_plex_box_faces 2,2 -dm_plex_separate_marker -bd 1,3 \
187: -dm_plex_transform_extrude_thickness 0.5
189: test:
190: suffix: tri_adapt_1
191: requires: triangle
192: nsize: 2
193: args: -dm_plex_box_faces 2,2 -bd 1 \
194: -dm_plex_transform_extrude_thickness 0.5 -petscpartitioner_type simple
196: test:
197: suffix: quad_adapt_0
198: args: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -dm_plex_separate_marker -bd 1,3 \
199: -dm_plex_transform_extrude_thickness 0.5
201: test:
202: suffix: quad_adapt_1
203: nsize: 2
204: args: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -bd 1 \
205: -dm_plex_transform_extrude_thickness 0.5 -petscpartitioner_type simple
207: test:
208: suffix: quad_adapt_1_forest
209: requires: p4est
210: nsize: 2
211: args: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -bd 1 \
212: -dm_plex_transform_extrude_thickness 0.5 -petscpartitioner_type simple \
213: -test_forest 1
215: test:
216: suffix: tet_adapt_0
217: requires: ctetgen
218: args: -dm_plex_dim 3 -dm_plex_box_faces 2,2,2 -dm_plex_separate_marker -bd 1,3 \
219: -dm_plex_transform_extrude_thickness 0.5
221: test:
222: suffix: tet_adapt_1
223: requires: ctetgen
224: nsize: 2
225: args: -dm_plex_dim 3 -dm_plex_box_faces 2,2,2 -bd 1 \
226: -dm_plex_transform_extrude_thickness 0.5 -petscpartitioner_type simple
228: test:
229: suffix: hex_adapt_0
230: args: -dm_plex_simplex 0 -dm_plex_dim 3 -dm_plex_box_faces 2,2,2 -dm_plex_separate_marker -bd 1,3 \
231: -dm_plex_transform_extrude_thickness 0.5
233: TEST*/