1 static const char help[] = "Tests for mesh extrusion";
2
3 #include <petscdmplex.h>
4 #include <petscdmforest.h>
5
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;
12
13 PETSC_EXTERN PetscErrorCode pyramidNormal(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *);
14
15 /* The pyramid apex is at (0.5, 0.5, -1) */
pyramidNormal(PetscInt dim,PetscReal time,const PetscReal x[],PetscInt r,PetscScalar u[],PetscCtx ctx)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;
20
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 }
25
ProcessOptions(MPI_Comm comm,AppCtx * options)26 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
27 {
28 PetscInt n = 64;
29 PetscBool flg;
30
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 }
42
CreateMesh(MPI_Comm comm,AppCtx * ctx,DM * dm)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 }
65
CreateAdaptLabel(DM dm,AppCtx * ctx,DMLabel * adaptLabel)66 static PetscErrorCode CreateAdaptLabel(DM dm, AppCtx *ctx, DMLabel *adaptLabel)
67 {
68 DMLabel label;
69 PetscInt b;
70
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;
82
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 }
93
main(int argc,char ** argv)94 int main(int argc, char **argv)
95 {
96 DM dm, dma;
97 DMLabel adaptLabel;
98 AppCtx ctx;
99
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) {
106 PetscCall(DMAdaptLabel(dm, adaptLabel, &dma));
107 } else {
108 PetscCall(DMExtrude(dm, 3, &dma));
109 }
110 PetscCall(PetscObjectSetName((PetscObject)dma, "Adapted Mesh"));
111 PetscCall(DMLabelDestroy(&adaptLabel));
112 PetscCall(DMDestroy(&dm));
113 PetscCall(DMViewFromOptions(dma, NULL, "-adapt_dm_view"));
114 PetscCall(DMDestroy(&dma));
115 PetscCall(PetscFinalize());
116 return 0;
117 }
118
119 /*TEST
120
121 testset:
122 args: -dm_view -adapt_dm_view -dm_plex_check_all
123
124 test:
125 suffix: seg_periodic_0
126 args: -dm_plex_dim 1 -dm_plex_box_faces 3 -dm_plex_transform_extrude_periodic -dm_plex_transform_extrude_use_tensor 0
127
128 test:
129 suffix: tri_tensor_0
130 requires: triangle
131 args: -dm_plex_transform_extrude_use_tensor {{0 1}separate output}
132
133 test:
134 suffix: quad_tensor_0
135 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_use_tensor {{0 1}separate output}
136
137 test:
138 suffix: quad_tensor_0_forest
139 requires: p4est
140 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_use_tensor {{0 1}separate output} -test_forest 1
141
142 test:
143 suffix: quad_normal_0
144 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal 0,1,1
145
146 test:
147 suffix: quad_normal_0_forest
148 requires: p4est
149 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal 0,1,1 -test_forest 1
150
151 test:
152 suffix: quad_normal_1
153 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal_function pyramidNormal
154
155 test:
156 suffix: quad_normal_1_forest
157 requires: p4est
158 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal_function pyramidNormal -test_forest 1
159
160 test:
161 suffix: quad_symmetric_0
162 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_symmetric
163
164 test:
165 suffix: quad_symmetric_0_forest
166 requires: p4est
167 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_symmetric -test_forest 1
168
169 test:
170 suffix: quad_label
171 args: -dm_plex_simplex 0 -dm_plex_transform_label_replica_inc {{0 100}separate output}
172
173 test:
174 suffix: quad_periodic_0
175 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_periodic -dm_plex_transform_extrude_use_tensor 0
176
177 testset:
178 args: -dm_adaptor cellrefiner -dm_plex_transform_type extrude \
179 -dm_view -adapt_dm_view -dm_plex_check_all
180
181 test:
182 suffix: tri_adapt_0
183 requires: triangle
184 args: -dm_plex_box_faces 2,2 -dm_plex_separate_marker -bd 1,3 \
185 -dm_plex_transform_extrude_thickness 0.5
186
187 test:
188 suffix: tri_adapt_1
189 requires: triangle
190 nsize: 2
191 args: -dm_plex_box_faces 2,2 -bd 1 \
192 -dm_plex_transform_extrude_thickness 0.5 -petscpartitioner_type simple
193
194 test:
195 suffix: quad_adapt_0
196 args: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -dm_plex_separate_marker -bd 1,3 \
197 -dm_plex_transform_extrude_thickness 0.5
198
199 test:
200 suffix: quad_adapt_1
201 nsize: 2
202 args: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -bd 1 \
203 -dm_plex_transform_extrude_thickness 0.5 -petscpartitioner_type simple
204
205 test:
206 suffix: quad_adapt_1_forest
207 requires: p4est
208 nsize: 2
209 args: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -bd 1 \
210 -dm_plex_transform_extrude_thickness 0.5 -petscpartitioner_type simple \
211 -test_forest 1
212
213 test:
214 suffix: tet_adapt_0
215 requires: ctetgen
216 args: -dm_plex_dim 3 -dm_plex_box_faces 2,2,2 -dm_plex_separate_marker -bd 1,3 \
217 -dm_plex_transform_extrude_thickness 0.5
218
219 test:
220 suffix: tet_adapt_1
221 requires: ctetgen
222 nsize: 2
223 args: -dm_plex_dim 3 -dm_plex_box_faces 2,2,2 -bd 1 \
224 -dm_plex_transform_extrude_thickness 0.5 -petscpartitioner_type simple
225
226 test:
227 suffix: hex_adapt_0
228 args: -dm_plex_simplex 0 -dm_plex_dim 3 -dm_plex_box_faces 2,2,2 -dm_plex_separate_marker -bd 1,3 \
229 -dm_plex_transform_extrude_thickness 0.5
230
231 TEST*/
232