1 static char help[] = "Tests for ephemeral meshes.\n"; 2 3 #include <petscdmplex.h> 4 #include <petscdmplextransform.h> 5 #include <petscdmlabelephemeral.h> 6 7 /* 8 Use 9 10 -ref_dm_view -eph_dm_view 11 12 to view the concrete and ephemeral meshes from the first transformation, and 13 14 -ref_dm_sec_view -eph_dm_sec_view 15 16 for the second. 17 */ 18 19 // Should remove when I have an API for everything 20 #include <petsc/private/dmplextransformimpl.h> 21 22 typedef struct { 23 DMLabel active; /* Label for transform */ 24 PetscBool second; /* Flag to execute a second transformation */ 25 PetscBool concrete; /* Flag to use the concrete mesh for the second transformation */ 26 } AppCtx; 27 28 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 29 { 30 PetscInt cells[1024], Nc = 1024; 31 PetscBool flg; 32 33 PetscFunctionBeginUser; 34 options->active = NULL; 35 options->second = PETSC_FALSE; 36 options->concrete = PETSC_FALSE; 37 38 PetscOptionsBegin(comm, "", "Ephemeral Meshing Options", "DMPLEX"); 39 PetscCall(PetscOptionsIntArray("-cells", "Cells to mark for transformation", "ex57.c", cells, &Nc, &flg)); 40 if (flg) { 41 PetscCall(DMLabelCreate(comm, "active", &options->active)); 42 for (PetscInt c = 0; c < Nc; ++c) PetscCall(DMLabelSetValue(options->active, cells[c], DM_ADAPT_REFINE)); 43 } 44 PetscCall(PetscOptionsBool("-second", "Use a second transformation", "ex57.c", options->second, &options->second, &flg)); 45 PetscCall(PetscOptionsBool("-concrete", "Use concrete mesh for the second transformation", "ex57.c", options->concrete, &options->concrete, &flg)); 46 PetscOptionsEnd(); 47 PetscFunctionReturn(PETSC_SUCCESS); 48 } 49 50 static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm) 51 { 52 PetscFunctionBeginUser; 53 PetscCall(DMCreate(comm, dm)); 54 PetscCall(DMSetType(*dm, DMPLEX)); 55 PetscCall(DMSetFromOptions(*dm)); 56 PetscCall(PetscObjectSetName((PetscObject)*dm, "Mesh")); 57 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 58 PetscFunctionReturn(PETSC_SUCCESS); 59 } 60 61 static PetscErrorCode CreateTransform(DM dm, DMLabel active, const char prefix[], DMPlexTransform *tr) 62 { 63 PetscFunctionBeginUser; 64 PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), tr)); 65 PetscCall(PetscObjectSetName((PetscObject)*tr, "Transform")); 66 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*tr, prefix)); 67 PetscCall(DMPlexTransformSetDM(*tr, dm)); 68 PetscCall(DMPlexTransformSetActive(*tr, active)); 69 PetscCall(DMPlexTransformSetFromOptions(*tr)); 70 PetscCall(DMPlexTransformSetUp(*tr)); 71 72 PetscCall(DMSetApplicationContext(dm, *tr)); 73 PetscCall(PetscObjectViewFromOptions((PetscObject)*tr, NULL, "-dm_plex_transform_view")); 74 PetscFunctionReturn(PETSC_SUCCESS); 75 } 76 77 static PetscErrorCode CreateEphemeralMesh(DMPlexTransform tr, DM *tdm) 78 { 79 PetscFunctionBegin; 80 PetscCall(DMPlexCreateEphemeral(tr, "eph_", tdm)); 81 PetscCall(PetscObjectSetName((PetscObject)*tdm, "Ephemeral Mesh")); 82 PetscFunctionReturn(PETSC_SUCCESS); 83 } 84 85 static PetscErrorCode CreateConcreteMesh(DMPlexTransform tr, DM *rdm) 86 { 87 DM cdm, codm, rcodm; 88 89 PetscFunctionBegin; 90 PetscCall(DMPlexTransformGetDM(tr, &cdm)); 91 PetscCall(DMPlexTransformApply(tr, cdm, rdm)); 92 PetscCall(DMSetCoarsenLevel(*rdm, cdm->leveldown)); 93 PetscCall(DMSetRefineLevel(*rdm, cdm->levelup + 1)); 94 PetscCall(DMCopyDisc(cdm, *rdm)); 95 PetscCall(DMGetCoordinateDM(cdm, &codm)); 96 PetscCall(DMGetCoordinateDM(*rdm, &rcodm)); 97 PetscCall(DMCopyDisc(codm, rcodm)); 98 PetscCall(DMPlexTransformCreateDiscLabels(tr, *rdm)); 99 PetscCall(DMSetCoarseDM(*rdm, cdm)); 100 PetscCall(DMPlexSetRegularRefinement(*rdm, PETSC_TRUE)); 101 if (rdm) { 102 ((DM_Plex *)(*rdm)->data)->printFEM = ((DM_Plex *)cdm->data)->printFEM; 103 ((DM_Plex *)(*rdm)->data)->printL2 = ((DM_Plex *)cdm->data)->printL2; 104 } 105 PetscCall(PetscObjectSetName((PetscObject)*rdm, "Concrete Mesh")); 106 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*rdm, "ref_")); 107 PetscFunctionReturn(PETSC_SUCCESS); 108 } 109 110 // dmA is concrete and dmB is ephemeral 111 static PetscErrorCode CompareMeshes(DM dmA, DM dmB, DM dm) 112 { 113 PetscInt dim, dimB, pStart, pEnd, pStartB, pEndB; 114 MPI_Comm comm; 115 116 PetscFunctionBegin; 117 PetscCall(PetscObjectGetComm((PetscObject)dmA, &comm)); 118 PetscCall(DMGetDimension(dmA, &dim)); 119 PetscCall(DMGetDimension(dmB, &dimB)); 120 PetscCheck(dim == dimB, comm, PETSC_ERR_ARG_INCOMP, "Dimension from dmA %" PetscInt_FMT " != %" PetscInt_FMT " from dmB", dim, dimB); 121 PetscCall(DMPlexGetChart(dmA, &pStart, &pEnd)); 122 PetscCall(DMPlexGetChart(dmB, &pStartB, &pEndB)); 123 PetscCheck(pStart == pStartB && pEnd == pEndB, comm, PETSC_ERR_ARG_INCOMP, "Chart from dmA (%" PetscInt_FMT ", %" PetscInt_FMT ") does not match (%" PetscInt_FMT ", %" PetscInt_FMT ") for dmB", pStart, pEnd, pStartB, pEndB); 124 for (PetscInt p = pStart; p < pEnd; ++p) { 125 const PetscInt *cone, *ornt, *coneB, *orntB; 126 PetscInt coneSize, coneSizeB; 127 128 PetscCall(DMPlexGetConeSize(dmA, p, &coneSize)); 129 PetscCall(DMPlexGetConeSize(dmB, p, &coneSizeB)); 130 PetscCheck(coneSize == coneSizeB, comm, PETSC_ERR_ARG_INCOMP, "Cone size for %" PetscInt_FMT " from dmA %" PetscInt_FMT " does not match %" PetscInt_FMT " for dmB", p, coneSize, coneSizeB); 131 PetscCall(DMPlexGetOrientedCone(dmA, p, &cone, &ornt)); 132 PetscCall(DMPlexGetOrientedCone(dmB, p, &coneB, &orntB)); 133 for (PetscInt c = 0; c < coneSize; ++c) { 134 PetscCheck(cone[c] == coneB[c] && ornt[c] == orntB[c], comm, PETSC_ERR_ARG_INCOMP, "Cone point %" PetscInt_FMT " for point %" PetscInt_FMT " from dmA (%" PetscInt_FMT ", %" PetscInt_FMT ") does not match (%" PetscInt_FMT ", %" PetscInt_FMT ") for dmB", c, p, cone[c], ornt[c], coneB[c], orntB[c]); 135 } 136 PetscCall(DMPlexRestoreOrientedCone(dmA, p, &cone, &ornt)); 137 PetscCall(DMPlexRestoreOrientedCone(dmB, p, &coneB, &orntB)); 138 } 139 PetscFunctionReturn(PETSC_SUCCESS); 140 } 141 142 int main(int argc, char *argv[]) 143 { 144 DM dm, tdm, rdm; 145 DMPlexTransform tr; 146 AppCtx user; 147 MPI_Comm comm; 148 149 PetscFunctionBeginUser; 150 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 151 comm = PETSC_COMM_WORLD; 152 PetscCall(ProcessOptions(comm, &user)); 153 PetscCall(CreateMesh(comm, &dm)); 154 PetscCall(CreateTransform(dm, user.active, "first_", &tr)); 155 PetscCall(CreateEphemeralMesh(tr, &tdm)); 156 PetscCall(CreateConcreteMesh(tr, &rdm)); 157 if (user.second) { 158 DMPlexTransform tr2; 159 DM tdm2, rdm2; 160 161 PetscCall(DMViewFromOptions(rdm, NULL, "-dm_sec_view")); 162 PetscCall(DMViewFromOptions(tdm, NULL, "-dm_sec_view")); 163 if (user.concrete) PetscCall(CreateTransform(rdm, user.active, "second_", &tr2)); 164 else PetscCall(CreateTransform(tdm, user.active, "second_", &tr2)); 165 PetscCall(CreateEphemeralMesh(tr2, &tdm2)); 166 PetscCall(CreateConcreteMesh(tr2, &rdm2)); 167 PetscCall(DMDestroy(&tdm)); 168 PetscCall(DMDestroy(&rdm)); 169 PetscCall(DMPlexTransformDestroy(&tr2)); 170 tdm = tdm2; 171 rdm = rdm2; 172 } 173 PetscCall(DMViewFromOptions(tdm, NULL, "-dm_view")); 174 PetscCall(DMViewFromOptions(rdm, NULL, "-dm_view")); 175 PetscCall(CompareMeshes(rdm, tdm, dm)); 176 PetscCall(DMPlexTransformDestroy(&tr)); 177 PetscCall(DMDestroy(&tdm)); 178 PetscCall(DMDestroy(&rdm)); 179 PetscCall(DMDestroy(&dm)); 180 PetscCall(DMLabelDestroy(&user.active)); 181 PetscCall(PetscFinalize()); 182 return 0; 183 } 184 185 /*TEST 186 187 # Tests for regular refinement of whole meshes 188 test: 189 suffix: tri 190 requires: triangle 191 args: -first_dm_plex_transform_view ::ascii_info_detail 192 193 test: 194 suffix: quad 195 args: -dm_plex_simplex 0 196 output_file: output/empty.out 197 198 # Here I am checking that the 'marker' label is correct for the ephemeral mesh 199 test: 200 suffix: tet 201 requires: ctetgen 202 args: -dm_plex_dim 3 -ref_dm_view -eph_dm_view -eph_dm_plex_view_labels marker 203 204 test: 205 suffix: hex 206 args: -dm_plex_dim 3 -dm_plex_simplex 0 207 output_file: output/empty.out 208 209 # Tests for filter patches 210 testset: 211 args: -first_dm_plex_transform_type transform_filter -ref_dm_view 212 213 test: 214 suffix: tri_patch 215 requires: triangle 216 args: -cells 0,1,2,4 217 218 test: 219 suffix: quad_patch 220 args: -dm_plex_simplex 0 -dm_plex_box_faces 3,3 -cells 0,1,3,4 221 222 # Tests for refined filter patches 223 testset: 224 args: -first_dm_plex_transform_type transform_filter -ref_dm_view -eph_dm_view -second 225 226 test: 227 suffix: tri_patch_ref 228 requires: triangle 229 args: -cells 0,1,2,4 230 231 test: 232 suffix: tri_patch_ref_concrete 233 requires: triangle 234 args: -cells 0,1,2,4 -concrete -first_dm_plex_transform_view ::ascii_info_detail 235 236 # Tests for boundary layer refinement 237 test: 238 suffix: quad_bl 239 args: -dm_plex_simplex 0 -dm_plex_dim 1 -dm_plex_box_faces 5 -dm_extrude 2 -cells 0,2,4,6,8 \ 240 -first_dm_plex_transform_type refine_boundary_layer -first_dm_plex_transform_bl_splits 4 \ 241 -ref_dm_view 242 243 # Tests for extrusion 244 test: 245 suffix: sphere_extruded 246 args: -dm_plex_shape sphere \ 247 -first_dm_plex_transform_type extrude \ 248 -first_dm_plex_transform_extrude_layers 3 \ 249 -first_dm_plex_transform_extrude_use_tensor 0 250 output_file: output/empty.out 251 252 TEST*/ 253