1 static const char help[] = "Tests for mesh extrusion"; 2 3 #include <petscdmplex.h> 4 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; 10 11 PETSC_EXTERN PetscErrorCode pyramidNormal(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar[], void *); 12 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; 18 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 } 23 24 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 25 { 26 PetscInt n = 64; 27 PetscBool flg; 28 PetscErrorCode ierr; 29 30 PetscFunctionBeginUser; 31 PetscCall(PetscStrcpy(options->bdLabel, "marker")); 32 ierr = PetscOptionsBegin(comm, "", "Parallel Mesh Adaptation Options", "DMPLEX");PetscCall(ierr); 33 PetscCall(PetscOptionsString("-label", "The boundary label name", "ex44.c", options->bdLabel, options->bdLabel, sizeof(options->bdLabel), NULL)); 34 PetscCall(PetscOptionsIntArray("-bd", "The boundaries to be extruded", "ex44.c", options->bd, &n, &flg)); 35 options->Nbd = flg ? n : 0; 36 ierr = PetscOptionsEnd();PetscCall(ierr); 37 PetscFunctionReturn(0); 38 } 39 40 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *ctx, DM *dm) 41 { 42 PetscFunctionBegin; 43 PetscCall(DMCreate(comm, dm)); 44 PetscCall(DMSetType(*dm, DMPLEX)); 45 PetscCall(DMSetFromOptions(*dm)); 46 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 47 PetscFunctionReturn(0); 48 } 49 50 static PetscErrorCode CreateAdaptLabel(DM dm, AppCtx *ctx, DMLabel *adaptLabel) 51 { 52 DMLabel label; 53 PetscInt b; 54 55 PetscFunctionBegin; 56 if (!ctx->Nbd) {*adaptLabel = NULL; PetscFunctionReturn(0);} 57 PetscCall(DMGetLabel(dm, ctx->bdLabel, &label)); 58 PetscCall(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; 63 64 PetscCall(DMLabelGetStratumIS(label, ctx->bd[b], &bdIS)); 65 if (!bdIS) continue; 66 PetscCall(ISGetLocalSize(bdIS, &n)); 67 PetscCall(ISGetIndices(bdIS, &points)); 68 for (i = 0; i < n; ++i) {PetscCall(DMLabelSetValue(*adaptLabel, points[i], DM_ADAPT_REFINE));} 69 PetscCall(ISRestoreIndices(bdIS, &points)); 70 PetscCall(ISDestroy(&bdIS)); 71 } 72 PetscFunctionReturn(0); 73 } 74 75 int main(int argc, char **argv) 76 { 77 DM dm, dma; 78 DMLabel adaptLabel; 79 AppCtx ctx; 80 81 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 82 PetscCall(ProcessOptions(PETSC_COMM_WORLD, &ctx)); 83 PetscCall(CreateMesh(PETSC_COMM_WORLD, &ctx, &dm)); 84 PetscCall(CreateAdaptLabel(dm, &ctx, &adaptLabel)); 85 if (adaptLabel) {PetscCall(DMAdaptLabel(dm, adaptLabel, &dma));} 86 else {PetscCall(DMExtrude(dm, 3, &dma));} 87 PetscCall(PetscObjectSetName((PetscObject) dma, "Adapted Mesh")); 88 PetscCall(DMLabelDestroy(&adaptLabel)); 89 PetscCall(DMDestroy(&dm)); 90 PetscCall(DMViewFromOptions(dma, NULL, "-adapt_dm_view")); 91 PetscCall(DMDestroy(&dma)); 92 PetscCall(PetscFinalize()); 93 return 0; 94 } 95 96 /*TEST 97 98 test: 99 suffix: tri_tensor_0 100 requires: triangle 101 args: -dm_plex_transform_extrude_use_tensor {{0 1}separate output} \ 102 -dm_view -adapt_dm_view -dm_plex_check_all 103 104 test: 105 suffix: quad_tensor_0 106 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_use_tensor {{0 1}separate output} \ 107 -dm_view -adapt_dm_view -dm_plex_check_all 108 109 test: 110 suffix: quad_normal_0 111 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal 0,1,1 \ 112 -dm_view -adapt_dm_view -dm_plex_check_all 113 114 test: 115 suffix: quad_normal_1 116 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_normal_function pyramidNormal \ 117 -dm_view -adapt_dm_view -dm_plex_check_all 118 119 test: 120 suffix: quad_symmetric_0 121 args: -dm_plex_simplex 0 -dm_plex_transform_extrude_symmetric \ 122 -dm_view -adapt_dm_view -dm_plex_check_all 123 124 testset: 125 args: -dm_adaptor cellrefiner -dm_plex_transform_type extrude \ 126 -dm_view -adapt_dm_view 127 128 test: 129 suffix: quad_adapt_0 130 args: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -dm_plex_separate_marker -bd 1,3 131 132 TEST*/ 133