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