Actual source code: ex41.c
1: static const char help[] = "Tests for adaptive refinement";
3: #include <petscdmplex.h>
5: typedef struct {
6: PetscBool metric; /* Flag to use metric adaptation, instead of tagging */
7: PetscInt *refcell; /* A cell to be refined on each process */
8: } AppCtx;
10: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
11: {
12: PetscMPIInt size;
13: PetscInt n;
16: options->metric = PETSC_FALSE;
17: MPI_Comm_size(comm, &size);
18: PetscCalloc1(size, &options->refcell);
19: n = size;
21: PetscOptionsBegin(comm, "", "Parallel Mesh Adaptation Options", "DMPLEX");
22: PetscOptionsBool("-metric", "Flag for metric refinement", "ex41.c", options->metric, &options->metric, NULL);
23: PetscOptionsIntArray("-refcell", "The cell to be refined", "ex41.c", options->refcell, &n, NULL);
25: PetscOptionsEnd();
26: return 0;
27: }
29: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *ctx, DM *dm)
30: {
31: DMCreate(comm, dm);
32: DMSetType(*dm, DMPLEX);
33: DMSetFromOptions(*dm);
34: DMViewFromOptions(*dm, NULL, "-dm_view");
35: return 0;
36: }
38: static PetscErrorCode CreateAdaptLabel(DM dm, AppCtx *ctx, DMLabel *adaptLabel)
39: {
40: PetscMPIInt rank;
42: MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank);
43: DMLabelCreate(PETSC_COMM_SELF, "Adaptation Label", adaptLabel);
44: if (ctx->refcell[rank] >= 0) DMLabelSetValue(*adaptLabel, ctx->refcell[rank], DM_ADAPT_REFINE);
45: return 0;
46: }
48: int main(int argc, char **argv)
49: {
50: DM dm, dma;
51: DMLabel adaptLabel;
52: AppCtx ctx;
55: PetscInitialize(&argc, &argv, NULL, help);
56: ProcessOptions(PETSC_COMM_WORLD, &ctx);
57: CreateMesh(PETSC_COMM_WORLD, &ctx, &dm);
58: CreateAdaptLabel(dm, &ctx, &adaptLabel);
59: DMAdaptLabel(dm, adaptLabel, &dma);
60: PetscObjectSetName((PetscObject)dma, "Adapted Mesh");
61: DMLabelDestroy(&adaptLabel);
62: DMDestroy(&dm);
63: DMViewFromOptions(dma, NULL, "-adapt_dm_view");
64: DMDestroy(&dma);
65: PetscFree(ctx.refcell);
66: PetscFinalize();
67: return 0;
68: }
70: /*TEST
72: testset:
73: args: -dm_adaptor cellrefiner -dm_plex_transform_type refine_sbr
75: test:
76: suffix: 0
77: requires: triangle
78: args: -dm_view -adapt_dm_view
80: test:
81: suffix: 1
82: requires: triangle
83: args: -dm_coord_space 0 -refcell 2 -dm_view ::ascii_info_detail -adapt_dm_view ::ascii_info_detail
85: test:
86: suffix: 2
87: requires: triangle
88: nsize: 2
89: args: -refcell 2,-1 -petscpartitioner_type simple -dm_view -adapt_dm_view
91: TEST*/