Actual source code: ex44.c

  1: static const char help[] = "Tests for mesh extrusion";

  3: #include <petscdmplex.h>

  5: typedef struct {
  6:   char     bdLabel[PETSC_MAX_PATH_LEN]; /* The boundary label name */
  7:   PetscInt Nbd;                         /* The number of boundary markers to extrude, 0 for all */
  8:   PetscInt bd[64];                      /* The boundary markers to be extruded */
  9: } AppCtx;

 11: PETSC_EXTERN PetscErrorCode pyramidNormal(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *);

 13: /* The pyramid apex is at (0.5, 0.5, -1) */
 14: PetscErrorCode pyramidNormal(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt r, PetscScalar u[], void *ctx)
 15: {
 16:   PetscReal apex[3] = {0.5, 0.5, -1.0};
 17:   PetscInt  d;

 19:   for (d = 0; d < dim; ++d) u[d] = x[d] - apex[d];
 20:   for (d = dim; d < 3; ++d) u[d] = 0.0 - apex[d];
 21:   return 0;
 22: }

 24: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 25: {
 26:   PetscInt  n = 64;
 27:   PetscBool flg;

 30:   PetscStrcpy(options->bdLabel, "marker");
 31:   PetscOptionsBegin(comm, "", "Parallel Mesh Adaptation Options", "DMPLEX");
 32:   PetscOptionsString("-label", "The boundary label name", "ex44.c", options->bdLabel, options->bdLabel, sizeof(options->bdLabel), NULL);
 33:   PetscOptionsIntArray("-bd", "The boundaries to be extruded", "ex44.c", options->bd, &n, &flg);
 34:   options->Nbd = flg ? n : 0;
 35:   PetscOptionsEnd();
 36:   return 0;
 37: }

 39: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *ctx, DM *dm)
 40: {
 41:   DMCreate(comm, dm);
 42:   DMSetType(*dm, DMPLEX);
 43:   DMSetFromOptions(*dm);
 44:   DMViewFromOptions(*dm, NULL, "-dm_view");
 45:   return 0;
 46: }

 48: static PetscErrorCode CreateAdaptLabel(DM dm, AppCtx *ctx, DMLabel *adaptLabel)
 49: {
 50:   DMLabel  label;
 51:   PetscInt b;

 53:   if (!ctx->Nbd) {
 54:     *adaptLabel = NULL;
 55:     return 0;
 56:   }
 57:   DMGetLabel(dm, ctx->bdLabel, &label);
 58:   DMLabelCreate(PETSC_COMM_SELF, "Adaptation Label", adaptLabel);
 59:   for (b = 0; b < ctx->Nbd; ++b) {
 60:     IS              bdIS;
 61:     const PetscInt *points;
 62:     PetscInt        n, i;

 64:     DMLabelGetStratumIS(label, ctx->bd[b], &bdIS);
 65:     if (!bdIS) continue;
 66:     ISGetLocalSize(bdIS, &n);
 67:     ISGetIndices(bdIS, &points);
 68:     for (i = 0; i < n; ++i) DMLabelSetValue(*adaptLabel, points[i], DM_ADAPT_REFINE);
 69:     ISRestoreIndices(bdIS, &points);
 70:     ISDestroy(&bdIS);
 71:   }
 72:   return 0;
 73: }

 75: int main(int argc, char **argv)
 76: {
 77:   DM      dm, dma;
 78:   DMLabel adaptLabel;
 79:   AppCtx  ctx;

 82:   PetscInitialize(&argc, &argv, NULL, help);
 83:   ProcessOptions(PETSC_COMM_WORLD, &ctx);
 84:   CreateMesh(PETSC_COMM_WORLD, &ctx, &dm);
 85:   CreateAdaptLabel(dm, &ctx, &adaptLabel);
 86:   if (adaptLabel) {
 87:     DMAdaptLabel(dm, adaptLabel, &dma);
 88:   } else {
 89:     DMExtrude(dm, 3, &dma);
 90:   }
 91:   PetscObjectSetName((PetscObject)dma, "Adapted Mesh");
 92:   DMLabelDestroy(&adaptLabel);
 93:   DMDestroy(&dm);
 94:   DMViewFromOptions(dma, NULL, "-adapt_dm_view");
 95:   DMDestroy(&dma);
 96:   PetscFinalize();
 97:   return 0;
 98: }

100: /*TEST

102:   test:
103:     suffix: tri_tensor_0
104:     requires: triangle
105:     args: -dm_plex_transform_extrude_use_tensor {{0 1}separate output} \
106:           -dm_view -adapt_dm_view -dm_plex_check_all

108:   test:
109:     suffix: quad_tensor_0
110:     args: -dm_plex_simplex 0 -dm_plex_transform_extrude_use_tensor {{0 1}separate output} \
111:           -dm_view -adapt_dm_view -dm_plex_check_all

113:   test:
114:     suffix: quad_normal_0
115:     args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal 0,1,1 \
116:           -dm_view -adapt_dm_view -dm_plex_check_all

118:   test:
119:     suffix: quad_normal_1
120:     args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal_function pyramidNormal \
121:           -dm_view -adapt_dm_view -dm_plex_check_all

123:   test:
124:     suffix: quad_symmetric_0
125:     args: -dm_plex_simplex 0 -dm_plex_transform_extrude_symmetric \
126:           -dm_view -adapt_dm_view -dm_plex_check_all

128:   testset:
129:     args: -dm_adaptor cellrefiner -dm_plex_transform_type extrude \
130:           -dm_view -adapt_dm_view

132:     test:
133:       suffix: quad_adapt_0
134:       args: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -dm_plex_separate_marker -bd 1,3 \
135:             -dm_plex_transform_extrude_thickness 0.5

137:     test:
138:       suffix: tet_adapt_0
139:       requires: ctetgen
140:       args: -dm_plex_dim 3 -dm_plex_box_faces 2,2,2 -dm_plex_separate_marker -bd 1,3 \
141:             -dm_plex_transform_extrude_thickness 0.5

143:     test:
144:       suffix: hex_adapt_0
145:       args: -dm_plex_simplex 0 -dm_plex_dim 3 -dm_plex_box_faces 2,2,2 -dm_plex_separate_marker -bd 1,3 \
146:             -dm_plex_transform_extrude_thickness 0.5

148: TEST*/