1 static char help[] = "Create a mesh, refine and coarsen simultaneously, and transfer a field\n\n";
2
3 #include <petscds.h>
4 #include <petscdmplex.h>
5 #include <petscdmforest.h>
6 #include <petscoptions.h>
7
AddIdentityLabel(DM dm)8 static PetscErrorCode AddIdentityLabel(DM dm)
9 {
10 PetscInt pStart, pEnd, p;
11
12 PetscFunctionBegin;
13 PetscCall(DMCreateLabel(dm, "identity"));
14 PetscCall(DMPlexGetChart(dm, &pStart, &pEnd));
15 for (p = pStart; p < pEnd; p++) PetscCall(DMSetLabelValue(dm, "identity", p, p));
16 PetscFunctionReturn(PETSC_SUCCESS);
17 }
18
CreateAdaptivityLabel(DM forest,DMLabel * adaptLabel)19 static PetscErrorCode CreateAdaptivityLabel(DM forest, DMLabel *adaptLabel)
20 {
21 DMLabel identLabel;
22 PetscInt cStart, cEnd, c;
23
24 PetscFunctionBegin;
25 PetscCall(DMLabelCreate(PETSC_COMM_SELF, "adapt", adaptLabel));
26 PetscCall(DMLabelSetDefaultValue(*adaptLabel, DM_ADAPT_COARSEN));
27 PetscCall(DMGetLabel(forest, "identity", &identLabel));
28 PetscCall(DMForestGetCellChart(forest, &cStart, &cEnd));
29 for (c = cStart; c < cEnd; c++) {
30 PetscInt basePoint;
31
32 PetscCall(DMLabelGetValue(identLabel, c, &basePoint));
33 if (!basePoint) PetscCall(DMLabelSetValue(*adaptLabel, c, DM_ADAPT_REFINE));
34 }
35 PetscFunctionReturn(PETSC_SUCCESS);
36 }
37
LinearFunction(PetscInt dim,PetscReal time,const PetscReal x[],PetscInt Nf,PetscScalar u[],PetscCtx ctx)38 static PetscErrorCode LinearFunction(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], PetscCtx ctx)
39 {
40 PetscFunctionBeginUser;
41 u[0] = (x[0] * 2.0 + 1.) + (x[1] * 20.0 + 10.) + ((dim == 3) ? (x[2] * 200.0 + 100.) : 0.);
42 PetscFunctionReturn(PETSC_SUCCESS);
43 }
44
MultiaffineFunction(PetscInt dim,PetscReal time,const PetscReal x[],PetscInt Nf,PetscScalar u[],PetscCtx ctx)45 static PetscErrorCode MultiaffineFunction(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], PetscCtx ctx)
46 {
47 PetscFunctionBeginUser;
48 u[0] = (x[0] * 1.0 + 2.0) * (x[1] * 3.0 - 4.0) * ((dim == 3) ? (x[2] * 5.0 + 6.0) : 1.);
49 PetscFunctionReturn(PETSC_SUCCESS);
50 }
51
CoordsFunction(PetscInt dim,PetscReal time,const PetscReal x[],PetscInt Nf,PetscScalar u[],PetscCtx ctx)52 static PetscErrorCode CoordsFunction(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], PetscCtx ctx)
53 {
54 PetscInt f;
55
56 PetscFunctionBeginUser;
57 for (f = 0; f < Nf; f++) u[f] = x[f];
58 PetscFunctionReturn(PETSC_SUCCESS);
59 }
60
61 typedef struct _bc_func_ctx {
62 PetscErrorCode (*func)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *);
63 PetscInt dim;
64 PetscInt Nf;
65 void *ctx;
66 } bc_func_ctx;
67
bc_func_fv(PetscReal time,const PetscReal * c,const PetscReal * n,const PetscScalar * xI,PetscScalar * xG,PetscCtx ctx)68 static PetscErrorCode bc_func_fv(PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *xI, PetscScalar *xG, PetscCtx ctx)
69 {
70 bc_func_ctx *bcCtx;
71
72 PetscFunctionBegin;
73 bcCtx = (bc_func_ctx *)ctx;
74 PetscCall(bcCtx->func(bcCtx->dim, time, c, bcCtx->Nf, xG, bcCtx->ctx));
75 PetscFunctionReturn(PETSC_SUCCESS);
76 }
77
IdentifyBadPoints(DM dm,Vec vec,PetscReal tol)78 static PetscErrorCode IdentifyBadPoints(DM dm, Vec vec, PetscReal tol)
79 {
80 DM dmplex;
81 PetscInt p, pStart, pEnd, maxDof;
82 Vec vecLocal;
83 DMLabel depthLabel;
84 PetscSection section;
85
86 PetscFunctionBegin;
87 PetscCall(DMCreateLocalVector(dm, &vecLocal));
88 PetscCall(DMGlobalToLocalBegin(dm, vec, INSERT_VALUES, vecLocal));
89 PetscCall(DMGlobalToLocalEnd(dm, vec, INSERT_VALUES, vecLocal));
90 PetscCall(DMConvert(dm, DMPLEX, &dmplex));
91 PetscCall(DMPlexGetChart(dmplex, &pStart, &pEnd));
92 PetscCall(DMPlexGetDepthLabel(dmplex, &depthLabel));
93 PetscCall(DMGetLocalSection(dmplex, §ion));
94 PetscCall(PetscSectionGetMaxDof(section, &maxDof));
95 for (p = pStart; p < pEnd; p++) {
96 PetscInt s, c, cSize, parent, childID, numChildren;
97 PetscInt cl, closureSize, *closure = NULL;
98 PetscScalar *values = NULL;
99 PetscBool bad = PETSC_FALSE;
100
101 PetscCall(VecGetValuesSection(vecLocal, section, p, &values));
102 PetscCall(PetscSectionGetDof(section, p, &cSize));
103 for (c = 0; c < cSize; c++) {
104 PetscReal absDiff = PetscAbsScalar(values[c]);
105 if (absDiff > tol) {
106 bad = PETSC_TRUE;
107 break;
108 }
109 }
110 if (!bad) continue;
111 PetscCall(PetscPrintf(PETSC_COMM_SELF, "Bad point %" PetscInt_FMT "\n", p));
112 PetscCall(DMLabelGetValue(depthLabel, p, &s));
113 PetscCall(PetscPrintf(PETSC_COMM_SELF, " Depth %" PetscInt_FMT "\n", s));
114 PetscCall(DMPlexGetTransitiveClosure(dmplex, p, PETSC_TRUE, &closureSize, &closure));
115 for (cl = 0; cl < closureSize; cl++) {
116 PetscInt cp = closure[2 * cl];
117 PetscCall(DMPlexGetTreeParent(dmplex, cp, &parent, &childID));
118 if (parent != cp) PetscCall(PetscPrintf(PETSC_COMM_SELF, " Closure point %" PetscInt_FMT " (%" PetscInt_FMT ") child of %" PetscInt_FMT " (ID %" PetscInt_FMT ")\n", cl, cp, parent, childID));
119 PetscCall(DMPlexGetTreeChildren(dmplex, cp, &numChildren, NULL));
120 if (numChildren) PetscCall(PetscPrintf(PETSC_COMM_SELF, " Closure point %" PetscInt_FMT " (%" PetscInt_FMT ") is parent\n", cl, cp));
121 }
122 PetscCall(DMPlexRestoreTransitiveClosure(dmplex, p, PETSC_TRUE, &closureSize, &closure));
123 for (c = 0; c < cSize; c++) {
124 PetscReal absDiff = PetscAbsScalar(values[c]);
125 if (absDiff > tol) PetscCall(PetscPrintf(PETSC_COMM_SELF, " Bad dof %" PetscInt_FMT "\n", c));
126 }
127 }
128 PetscCall(DMDestroy(&dmplex));
129 PetscCall(VecDestroy(&vecLocal));
130 PetscFunctionReturn(PETSC_SUCCESS);
131 }
132
main(int argc,char ** argv)133 int main(int argc, char **argv)
134 {
135 MPI_Comm comm;
136 DM base, preForest, postForest;
137 PetscInt dim, Nf = 1;
138 PetscInt step, adaptSteps = 1;
139 PetscInt preCount, postCount;
140 Vec preVec, postVecTransfer, postVecExact;
141 PetscErrorCode (*funcs[1])(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *) = {MultiaffineFunction};
142 void *ctxs[1] = {NULL};
143 PetscReal diff, tol = PETSC_SMALL;
144 PetscBool linear = PETSC_FALSE;
145 PetscBool coords = PETSC_FALSE;
146 PetscBool useFV = PETSC_FALSE;
147 PetscBool conv = PETSC_FALSE;
148 PetscBool transfer_from_base[2] = {PETSC_TRUE, PETSC_FALSE};
149 PetscBool use_bcs = PETSC_TRUE;
150 bc_func_ctx bcCtx;
151 DMLabel adaptLabel;
152
153 PetscFunctionBeginUser;
154 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
155 comm = PETSC_COMM_WORLD;
156 PetscOptionsBegin(comm, "", "DMForestTransferVec() Test Options", "DMFOREST");
157 PetscCall(PetscOptionsBool("-linear", "Transfer a simple linear function", "ex2.c", linear, &linear, NULL));
158 PetscCall(PetscOptionsBool("-coords", "Transfer a simple coordinate function", "ex2.c", coords, &coords, NULL));
159 PetscCall(PetscOptionsBool("-use_fv", "Use a finite volume approximation", "ex2.c", useFV, &useFV, NULL));
160 PetscCall(PetscOptionsBool("-test_convert", "Test conversion to DMPLEX", NULL, conv, &conv, NULL));
161 PetscCall(PetscOptionsBool("-transfer_from_base", "Transfer a vector from base DM to DMForest", "ex2.c", transfer_from_base[0], &transfer_from_base[0], NULL));
162 transfer_from_base[1] = transfer_from_base[0];
163 PetscCall(PetscOptionsBool("-transfer_from_base_steps", "Transfer a vector from base DM to the latest DMForest after the adaptivity steps", "ex2.c", transfer_from_base[1], &transfer_from_base[1], NULL));
164 PetscCall(PetscOptionsBool("-use_bcs", "Use dirichlet boundary conditions", "ex2.c", use_bcs, &use_bcs, NULL));
165 PetscCall(PetscOptionsBoundedInt("-adapt_steps", "Number of adaptivity steps", "ex2.c", adaptSteps, &adaptSteps, NULL, 0));
166 PetscOptionsEnd();
167
168 tol = PetscMax(1.e-10, tol); /* XXX fix for quadruple precision -> why do I need to do this? */
169
170 /* the base mesh */
171 PetscCall(DMCreate(comm, &base));
172 PetscCall(DMSetType(base, DMPLEX));
173 PetscCall(DMSetFromOptions(base));
174
175 PetscCall(AddIdentityLabel(base));
176 PetscCall(DMGetDimension(base, &dim));
177
178 if (linear) funcs[0] = LinearFunction;
179 if (coords) {
180 funcs[0] = CoordsFunction;
181 Nf = dim;
182 }
183
184 bcCtx.func = funcs[0];
185 bcCtx.dim = dim;
186 bcCtx.Nf = Nf;
187 bcCtx.ctx = NULL;
188
189 if (useFV) {
190 PetscFV fv;
191 PetscLimiter limiter;
192 DM baseFV;
193
194 PetscCall(DMPlexConstructGhostCells(base, NULL, NULL, &baseFV));
195 PetscCall(DMViewFromOptions(baseFV, NULL, "-fv_dm_view"));
196 PetscCall(DMDestroy(&base));
197 base = baseFV;
198 PetscCall(PetscFVCreate(PETSC_COMM_SELF, &fv));
199 PetscCall(PetscFVSetSpatialDimension(fv, dim));
200 PetscCall(PetscFVSetType(fv, PETSCFVLEASTSQUARES));
201 PetscCall(PetscFVSetNumComponents(fv, Nf));
202 PetscCall(PetscLimiterCreate(comm, &limiter));
203 PetscCall(PetscLimiterSetType(limiter, PETSCLIMITERNONE));
204 PetscCall(PetscFVSetLimiter(fv, limiter));
205 PetscCall(PetscLimiterDestroy(&limiter));
206 PetscCall(PetscFVSetFromOptions(fv));
207 PetscCall(DMSetField(base, 0, NULL, (PetscObject)fv));
208 PetscCall(PetscFVDestroy(&fv));
209 } else {
210 PetscFE fe;
211
212 PetscCall(PetscFECreateDefault(comm, dim, Nf, PETSC_FALSE, NULL, PETSC_DEFAULT, &fe));
213 PetscCall(DMSetField(base, 0, NULL, (PetscObject)fe));
214 PetscCall(PetscFEDestroy(&fe));
215 }
216 PetscCall(DMCreateDS(base));
217
218 if (use_bcs) {
219 PetscInt ids[] = {1, 2, 3, 4, 5, 6};
220 DMLabel label;
221
222 PetscCall(DMGetLabel(base, "marker", &label));
223 PetscCall(DMAddBoundary(base, DM_BC_ESSENTIAL, "bc", label, 2 * dim, ids, 0, 0, NULL, useFV ? (PetscVoidFn *)bc_func_fv : (PetscVoidFn *)funcs[0], NULL, useFV ? (void *)&bcCtx : NULL, NULL));
224 }
225 PetscCall(DMViewFromOptions(base, NULL, "-dm_base_view"));
226
227 /* the pre adaptivity forest */
228 PetscCall(DMCreate(comm, &preForest));
229 PetscCall(DMSetType(preForest, (dim == 2) ? DMP4EST : DMP8EST));
230 PetscCall(DMCopyDisc(base, preForest));
231 PetscCall(DMForestSetBaseDM(preForest, base));
232 PetscCall(DMForestSetMinimumRefinement(preForest, 0));
233 PetscCall(DMForestSetInitialRefinement(preForest, 1));
234 PetscCall(DMSetFromOptions(preForest));
235 PetscCall(DMSetUp(preForest));
236 PetscCall(DMViewFromOptions(preForest, NULL, "-dm_pre_view"));
237
238 /* the pre adaptivity field */
239 PetscCall(DMCreateGlobalVector(preForest, &preVec));
240 PetscCall(DMProjectFunction(preForest, 0., funcs, ctxs, INSERT_VALUES, preVec));
241 PetscCall(VecViewFromOptions(preVec, NULL, "-vec_pre_view"));
242
243 /* communicate between base and pre adaptivity forest */
244 if (transfer_from_base[0]) {
245 Vec baseVec, baseVecMapped;
246
247 PetscCall(DMGetGlobalVector(base, &baseVec));
248 PetscCall(DMProjectFunction(base, 0., funcs, ctxs, INSERT_VALUES, baseVec));
249 PetscCall(PetscObjectSetName((PetscObject)baseVec, "Function Base"));
250 PetscCall(VecViewFromOptions(baseVec, NULL, "-vec_base_view"));
251
252 PetscCall(DMGetGlobalVector(preForest, &baseVecMapped));
253 PetscCall(DMForestTransferVecFromBase(preForest, baseVec, baseVecMapped));
254 PetscCall(VecViewFromOptions(baseVecMapped, NULL, "-vec_map_base_view"));
255
256 /* compare */
257 PetscCall(VecAXPY(baseVecMapped, -1., preVec));
258 PetscCall(VecViewFromOptions(baseVecMapped, NULL, "-vec_map_diff_view"));
259 PetscCall(VecNorm(baseVecMapped, NORM_2, &diff));
260
261 /* output */
262 if (diff < tol) {
263 PetscCall(PetscPrintf(comm, "DMForestTransferVecFromBase() passes.\n"));
264 } else {
265 PetscCall(PetscPrintf(comm, "DMForestTransferVecFromBase() fails with error %g and tolerance %g\n", (double)diff, (double)tol));
266 }
267
268 PetscCall(DMRestoreGlobalVector(base, &baseVec));
269 PetscCall(DMRestoreGlobalVector(preForest, &baseVecMapped));
270 }
271
272 for (step = 0; step < adaptSteps; ++step) {
273 if (!transfer_from_base[1]) PetscCall(PetscObjectGetReference((PetscObject)preForest, &preCount));
274
275 /* adapt */
276 PetscCall(CreateAdaptivityLabel(preForest, &adaptLabel));
277 PetscCall(DMForestTemplate(preForest, comm, &postForest));
278 if (step) PetscCall(DMForestSetAdaptivityLabel(postForest, adaptLabel));
279 PetscCall(DMLabelDestroy(&adaptLabel));
280 PetscCall(DMSetUp(postForest));
281 PetscCall(DMViewFromOptions(postForest, NULL, "-dm_post_view"));
282
283 /* transfer */
284 PetscCall(DMCreateGlobalVector(postForest, &postVecTransfer));
285 PetscCall(DMForestTransferVec(preForest, preVec, postForest, postVecTransfer, PETSC_TRUE, 0.0));
286 PetscCall(VecViewFromOptions(postVecTransfer, NULL, "-vec_post_transfer_view"));
287
288 /* the exact post adaptivity field */
289 PetscCall(DMCreateGlobalVector(postForest, &postVecExact));
290 PetscCall(DMProjectFunction(postForest, 0., funcs, ctxs, INSERT_VALUES, postVecExact));
291 PetscCall(VecViewFromOptions(postVecExact, NULL, "-vec_post_exact_view"));
292
293 /* compare */
294 PetscCall(VecAXPY(postVecExact, -1., postVecTransfer));
295 PetscCall(VecViewFromOptions(postVecExact, NULL, "-vec_diff_view"));
296 PetscCall(VecNorm(postVecExact, NORM_2, &diff));
297
298 /* output */
299 if (diff < tol) {
300 PetscCall(PetscPrintf(comm, "DMForestTransferVec() passes.\n"));
301 } else {
302 PetscCall(PetscPrintf(comm, "DMForestTransferVec() fails with error %g and tolerance %g\n", (double)diff, (double)tol));
303 PetscCall(IdentifyBadPoints(postForest, postVecExact, tol));
304 }
305 PetscCall(VecDestroy(&postVecExact));
306
307 /* disconnect preForest from postForest if we don't test the transfer throughout the entire refinement process */
308 if (!transfer_from_base[1]) {
309 PetscCall(DMForestSetAdaptivityForest(postForest, NULL));
310 PetscCall(PetscObjectGetReference((PetscObject)preForest, &postCount));
311 PetscCheck(postCount == preCount, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Adaptation not memory neutral: reference count increase from %" PetscInt_FMT " to %" PetscInt_FMT, preCount, postCount);
312 }
313
314 if (conv) {
315 DM dmConv;
316
317 PetscCall(DMConvert(postForest, DMPLEX, &dmConv));
318 PetscCall(DMViewFromOptions(dmConv, NULL, "-dm_conv_view"));
319 PetscCall(DMPlexCheckCellShape(dmConv, PETSC_TRUE, PETSC_DETERMINE));
320 PetscCall(DMDestroy(&dmConv));
321 }
322
323 PetscCall(VecDestroy(&preVec));
324 PetscCall(DMDestroy(&preForest));
325
326 preVec = postVecTransfer;
327 preForest = postForest;
328 }
329
330 if (transfer_from_base[1]) {
331 Vec baseVec, baseVecMapped;
332
333 /* communicate between base and last adapted forest */
334 PetscCall(DMGetGlobalVector(base, &baseVec));
335 PetscCall(DMProjectFunction(base, 0., funcs, ctxs, INSERT_VALUES, baseVec));
336 PetscCall(PetscObjectSetName((PetscObject)baseVec, "Function Base"));
337 PetscCall(VecViewFromOptions(baseVec, NULL, "-vec_base_view"));
338
339 PetscCall(DMGetGlobalVector(preForest, &baseVecMapped));
340 PetscCall(DMForestTransferVecFromBase(preForest, baseVec, baseVecMapped));
341 PetscCall(VecViewFromOptions(baseVecMapped, NULL, "-vec_map_base_view"));
342
343 /* compare */
344 PetscCall(VecAXPY(baseVecMapped, -1., preVec));
345 PetscCall(VecViewFromOptions(baseVecMapped, NULL, "-vec_map_diff_view"));
346 PetscCall(VecNorm(baseVecMapped, NORM_2, &diff));
347
348 /* output */
349 if (diff < tol) {
350 PetscCall(PetscPrintf(comm, "DMForestTransferVecFromBase() passes.\n"));
351 } else {
352 PetscCall(PetscPrintf(comm, "DMForestTransferVecFromBase() fails with error %g and tolerance %g\n", (double)diff, (double)tol));
353 }
354
355 PetscCall(DMRestoreGlobalVector(base, &baseVec));
356 PetscCall(DMRestoreGlobalVector(preForest, &baseVecMapped));
357 }
358
359 /* cleanup */
360 PetscCall(VecDestroy(&preVec));
361 PetscCall(DMDestroy(&preForest));
362 PetscCall(DMDestroy(&base));
363 PetscCall(PetscFinalize());
364 return 0;
365 }
366
367 /*TEST
368 testset:
369 args: -dm_plex_simplex 0 -dm_plex_box_faces 3,3,3 -petscspace_type tensor
370
371 test:
372 output_file: output/ex2_2d.out
373 suffix: p4est_2d
374 args: -petscspace_degree 2
375 nsize: 3
376 requires: p4est !single
377
378 test:
379 output_file: output/ex2_2d.out
380 suffix: p4est_2d_deg4
381 args: -petscspace_degree 4
382 requires: p4est !single
383
384 test:
385 output_file: output/ex2_2d.out
386 suffix: p4est_2d_deg8
387 args: -petscspace_degree 8
388 requires: p4est !single
389
390 test:
391 output_file: output/ex2_steps2.out
392 suffix: p4est_2d_deg2_steps2
393 args: -petscspace_degree 2 -coords -adapt_steps 2
394 nsize: 3
395 requires: p4est !single
396
397 test:
398 output_file: output/ex2_steps3.out
399 suffix: p4est_2d_deg3_steps3
400 args: -petscspace_degree 3 -coords -adapt_steps 3 -petscdualspace_lagrange_node_type equispaced -petscdualspace_lagrange_node_endpoints 1
401 nsize: 3
402 requires: p4est !single
403
404 test:
405 output_file: output/ex2_steps3.out
406 suffix: p4est_2d_deg3_steps3_L2_periodic
407 args: -petscspace_degree 3 -petscdualspace_lagrange_continuity 0 -coords -adapt_steps 3 -dm_plex_box_bd periodic,periodic -use_bcs 0 -petscdualspace_lagrange_node_type equispaced
408 nsize: 3
409 requires: p4est !single
410
411 test:
412 output_file: output/ex2_steps3.out
413 suffix: p4est_3d_deg2_steps3_L2_periodic
414 args: -dm_plex_dim 3 -petscspace_degree 2 -petscdualspace_lagrange_continuity 0 -coords -adapt_steps 3 -dm_plex_box_bd periodic,periodic,periodic -use_bcs 0
415 nsize: 3
416 requires: p4est !single
417
418 test:
419 output_file: output/ex2_steps2.out
420 suffix: p4est_3d_deg2_steps2
421 args: -dm_plex_dim 3 -petscspace_degree 2 -coords -adapt_steps 2
422 nsize: 3
423 requires: p4est !single
424
425 test:
426 output_file: output/ex2_steps3.out
427 suffix: p4est_3d_deg3_steps3
428 args: -dm_plex_dim 3 -petscspace_degree 3 -coords -adapt_steps 3 -petscdualspace_lagrange_node_type equispaced -petscdualspace_lagrange_node_endpoints 1
429 nsize: 3
430 requires: p4est !single
431
432 test:
433 output_file: output/ex2_3d.out
434 suffix: p4est_3d
435 args: -dm_plex_dim 3 -petscspace_degree 1
436 nsize: 3
437 requires: p4est !single
438
439 test:
440 output_file: output/ex2_3d.out
441 suffix: p4est_3d_deg3
442 args: -dm_plex_dim 3 -petscspace_degree 3
443 nsize: 3
444 requires: p4est !single
445
446 test:
447 output_file: output/ex2_2d.out
448 suffix: p4est_2d_deg2_coords
449 args: -petscspace_degree 2 -coords
450 nsize: 3
451 requires: p4est !single
452
453 test:
454 output_file: output/ex2_3d.out
455 suffix: p4est_3d_deg2_coords
456 args: -dm_plex_dim 3 -petscspace_degree 2 -coords
457 nsize: 3
458 requires: p4est !single
459
460 test:
461 suffix: p4est_3d_nans
462 args: -dm_plex_dim 3 -dm_forest_partition_overlap 1 -test_convert -petscspace_degree 1
463 nsize: 2
464 requires: p4est !single
465
466 test:
467 TODO: not broken, but the 3D case below is broken, so I do not trust this one
468 output_file: output/ex2_steps2.out
469 suffix: p4est_2d_tfb_distributed_nc
470 args: -petscspace_degree 3 -dm_forest_maximum_refinement 2 -dm_p4est_refine_pattern hash -use_bcs 0 -coords -adapt_steps 2 -petscpartitioner_type shell -petscpartitioner_shell_random
471 nsize: 3
472 requires: p4est !single
473
474 test:
475 TODO: broken
476 output_file: output/ex2_steps2.out
477 suffix: p4est_3d_tfb_distributed_nc
478 args: -dm_plex_dim 3 -petscspace_degree 2 -dm_forest_maximum_refinement 2 -dm_p4est_refine_pattern hash -use_bcs 0 -coords -adapt_steps 2 -petscpartitioner_type shell -petscpartitioner_shell_random
479 nsize: 3
480 requires: p4est !single
481
482 testset:
483 args: -petscspace_type tensor -dm_coord_space 0 -dm_plex_transform_type refine_tobox
484
485 test:
486 TODO: broken
487 output_file: output/ex2_3d.out
488 suffix: p4est_3d_transfer_fails
489 args: -petscspace_degree 1 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -adapt_steps 1 -dm_forest_initial_refinement 1 -use_bcs 0 -dm_refine
490 requires: p4est !single
491
492 test:
493 TODO: broken
494 output_file: output/ex2_steps2_notfb.out
495 suffix: p4est_3d_transfer_fails_2
496 args: -petscspace_degree 1 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -adapt_steps 2 -dm_forest_initial_refinement 0 -transfer_from_base 0 -use_bcs 0 -dm_refine
497 requires: p4est !single
498
499 test:
500 output_file: output/ex2_steps2.out
501 suffix: p4est_3d_multi_transfer_s2t
502 args: -petscspace_degree 3 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -adapt_steps 2 -dm_forest_initial_refinement 1 -petscdualspace_lagrange_continuity 0 -use_bcs 0 -dm_refine 1
503 requires: p4est !single
504
505 test:
506 output_file: output/ex2_steps2.out
507 suffix: p4est_3d_coords_transfer_s2t
508 args: -petscspace_degree 3 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -adapt_steps 2 -dm_forest_initial_refinement 1 -petscdualspace_lagrange_continuity 0 -coords -use_bcs 0 -dm_refine 1
509 requires: p4est !single
510
511 testset:
512 args: -dm_plex_simplex 0 -dm_plex_box_faces 3,3,3
513
514 test:
515 output_file: output/ex2_2d_fv.out
516 suffix: p4est_2d_fv
517 args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1
518 nsize: 3
519 requires: p4est !single
520
521 test:
522 TODO: broken (codimension adjacency)
523 output_file: output/ex2_2d_fv.out
524 suffix: p4est_2d_fv_adjcodim
525 args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1 -dm_forest_adjacency_codimension 1
526 nsize: 2
527 requires: p4est !single
528
529 test:
530 TODO: broken (dimension adjacency)
531 output_file: output/ex2_2d_fv.out
532 suffix: p4est_2d_fv_adjdim
533 args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1 -dm_forest_adjacency_dimension 1
534 nsize: 2
535 requires: p4est !single
536
537 test:
538 output_file: output/ex2_2d_fv.out
539 suffix: p4est_2d_fv_zerocells
540 args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1
541 nsize: 10
542 requires: p4est !single
543
544 test:
545 output_file: output/ex2_3d_fv.out
546 suffix: p4est_3d_fv
547 args: -dm_plex_dim 3 -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1
548 nsize: 3
549 requires: p4est !single
550
551 TEST*/
552