Actual source code: ex41.c
1: static const char help[] = "Tests for adaptive refinement";
3: #include <petscdmplex.h>
4: #include <petscdmplextransform.h>
6: typedef struct {
7: PetscBool metric; /* Flag to use metric adaptation, instead of tagging */
8: PetscInt iter; /* Number of adaptation generations */
9: PetscInt *refcell; /* A cell to be refined on each process */
10: } AppCtx;
12: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
13: {
14: PetscMPIInt size;
15: PetscInt n;
17: PetscFunctionBeginUser;
18: options->metric = PETSC_FALSE;
19: options->iter = 1;
20: PetscCallMPI(MPI_Comm_size(comm, &size));
21: PetscCall(PetscCalloc1(size, &options->refcell));
22: n = size;
24: PetscOptionsBegin(comm, "", "Parallel Mesh Adaptation Options", "DMPLEX");
25: PetscCall(PetscOptionsBool("-metric", "Flag for metric refinement", "ex41.c", options->metric, &options->metric, NULL));
26: PetscCall(PetscOptionsInt("-adapt_iter", "Number of adaptation generations", "ex41.c", options->iter, &options->iter, NULL));
27: PetscCall(PetscOptionsIntArray("-refcell", "The cell to be refined", "ex41.c", options->refcell, &n, NULL));
28: if (n) PetscCheck(n == size, comm, PETSC_ERR_ARG_SIZ, "Only gave %" PetscInt_FMT " cells to refine, must give one for all %d processes", n, size);
29: PetscOptionsEnd();
30: PetscFunctionReturn(PETSC_SUCCESS);
31: }
33: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *ctx, DM *dm)
34: {
35: PetscFunctionBegin;
36: PetscCall(DMCreate(comm, dm));
37: PetscCall(DMSetType(*dm, DMPLEX));
38: PetscCall(DMSetFromOptions(*dm));
39: PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
40: PetscFunctionReturn(PETSC_SUCCESS);
41: }
43: static PetscErrorCode CreateAdaptLabel(DM dm, AppCtx *ctx, DMLabel *adaptLabel)
44: {
45: PetscMPIInt rank;
47: PetscFunctionBegin;
48: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
49: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Adaptation Label", adaptLabel));
50: if (ctx->refcell[rank] >= 0) PetscCall(DMLabelSetValue(*adaptLabel, ctx->refcell[rank], DM_ADAPT_REFINE));
51: PetscFunctionReturn(PETSC_SUCCESS);
52: }
54: static PetscErrorCode ConstructRefineTree(DM dm, DM odm)
55: {
56: DMPlexTransform tr;
57: PetscInt cStart, cEnd;
59: PetscFunctionBegin;
60: PetscCall(DMPlexGetTransform(dm, &tr));
61: if (!tr) PetscFunctionReturn(PETSC_SUCCESS);
62: PetscCall(DMPlexGetHeightStratum(odm, 0, &cStart, &cEnd));
63: for (PetscInt c = cStart; c < cEnd; ++c) {
64: DMPolytopeType ct;
65: DMPolytopeType *rct;
66: PetscInt *rsize, *rcone, *rornt;
67: PetscInt Nct, dim, pNew = 0;
69: PetscCall(PetscPrintf(PETSC_COMM_SELF, "Cell %" PetscInt_FMT " produced new cells", c));
70: PetscCall(DMPlexGetCellType(odm, c, &ct));
71: dim = DMPolytopeTypeGetDim(ct);
72: PetscCall(DMPlexTransformCellTransform(tr, ct, c, NULL, &Nct, &rct, &rsize, &rcone, &rornt));
73: for (PetscInt n = 0; n < Nct; ++n) {
74: if (DMPolytopeTypeGetDim(rct[n]) != dim) continue;
75: for (PetscInt r = 0; r < rsize[n]; ++r) {
76: PetscCall(DMPlexTransformGetTargetPoint(tr, ct, rct[n], c, r, &pNew));
77: PetscCall(PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT, pNew));
78: }
79: }
80: PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n"));
81: }
82: PetscFunctionReturn(PETSC_SUCCESS);
83: }
85: int main(int argc, char **argv)
86: {
87: DM dm, dma = NULL;
88: DMLabel adaptLabel;
89: AppCtx ctx;
91: PetscFunctionBeginUser;
92: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
93: PetscCall(ProcessOptions(PETSC_COMM_WORLD, &ctx));
94: PetscCall(CreateMesh(PETSC_COMM_WORLD, &ctx, &dm));
95: for (PetscInt i = 0; i < ctx.iter; ++i) {
96: PetscCall(CreateAdaptLabel(dm, &ctx, &adaptLabel));
97: PetscCall(DMAdaptLabel(dm, adaptLabel, &dma));
98: PetscCall(PetscObjectSetName((PetscObject)dma, "Adapted Mesh"));
99: PetscCall(DMLabelDestroy(&adaptLabel));
100: PetscCall(DMPlexCheck(dma));
101: PetscCall(DMViewFromOptions(dma, NULL, "-adapt_dm_view"));
102: if (i < ctx.iter - 1) {
103: PetscCall(DMDestroy(&dm));
104: dm = dma;
105: }
106: }
107: PetscCall(ConstructRefineTree(dma, dm));
108: PetscCall(DMDestroy(&dm));
109: PetscCall(DMDestroy(&dma));
110: PetscCall(PetscFree(ctx.refcell));
111: PetscCall(PetscFinalize());
112: return 0;
113: }
115: /*TEST
117: testset:
118: args: -dm_adaptor cellrefiner -dm_plex_transform_type refine_sbr
120: test:
121: suffix: 0
122: requires: triangle
123: args: -dm_view -adapt_dm_view
125: test:
126: suffix: 1
127: requires: triangle
128: args: -dm_coord_space 0 -refcell 2 -dm_view ::ascii_info_detail -adapt_dm_view ::ascii_info_detail
130: test:
131: suffix: 1_save
132: requires: triangle
133: args: -refcell 2 -dm_plex_save_transform -dm_view -adapt_dm_view
135: test:
136: suffix: 2
137: requires: triangle
138: nsize: 2
139: args: -refcell 2,-1 -petscpartitioner_type simple -dm_view -adapt_dm_view
141: # Refine one tetrahedron of two, the example of Fig. 6 of Plaza & Carey (2000), and validate
142: # the subdivision generator against the enumeration of Table 4 of the paper
143: test:
144: suffix: 3
145: args: -dm_plex_shape doublet -dm_plex_dim 3 -dm_plex_simplex 1 -refcell 0 \
146: -dm_plex_transform_sbr_validate -dm_view -adapt_dm_view
148: # Keep refining the first cell of each generation, verifying mesh validity each time
149: test:
150: suffix: 4
151: args: -dm_plex_shape doublet -dm_plex_dim 3 -dm_plex_simplex 1 -refcell 0 -adapt_iter 4 -adapt_dm_view
153: # One tetrahedron per process, exercising the conformity closure across the shared face
154: test:
155: suffix: 5
156: nsize: 2
157: args: -dm_plex_shape doublet -dm_plex_dim 3 -dm_plex_simplex 1 -refcell 0,-1 -adapt_iter 2 \
158: -petscpartitioner_type simple -dm_distribute -dm_view -adapt_dm_view
160: TEST*/