1 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 2 #include <petscdmplextransform.h> 3 4 /*@C 5 DMPlexExtrude - Extrude a volumetric mesh from the input surface mesh 6 7 Input Parameters: 8 + dm - The surface mesh 9 . layers - The number of extruded layers 10 . thickness - The total thickness of the extruded layers, or PETSC_DETERMINE 11 . tensor - Flag to create tensor produt cells 12 . symmetric - Flag to extrude symmetrically about the surface 13 . normal - Surface normal vector, or NULL 14 - thicknesses - Thickness of each layer, or NULL 15 16 Output Parameter: 17 . edm - The volumetric mesh 18 19 Notes: 20 Extrusion is implemented as a DMPlexTransform, so that new mesh points are produced from old mesh points. In the exmaple below, 21 we begin with an edge (v0, v3). It is extruded for two layers. The original vertex v0 produces two edges, e1 and e2, and three vertices, 22 v0, v2, and v2. Similarly, vertex v3 produces e3, e4, v3, v4, and v5. The original edge produces itself, e5 and e6, as well as face1 and 23 face2. The new mesh points are given the same labels as the original points which produced them. Thus, if v0 had a label value 1, then so 24 would v1, v2, e1 and e2. 25 26 $ v2----- e6 -----v5 27 $ | | 28 $ e2 face2 e4 29 $ | | 30 $ v1----- e5 -----v4 31 $ | | 32 $ e1 face1 e3 33 $ | | 34 $ v0--- original ----v3 35 36 Options Database: 37 + -dm_plex_transform_extrude_thickness <t> - The total thickness of extruded layers 38 . -dm_plex_transform_extrude_use_tensor <bool> - Use tensor cells when extruding 39 . -dm_plex_transform_extrude_symmetric <bool> - Extrude layers symmetrically about the surface 40 . -dm_plex_transform_extrude_normal <n0,...,nd> - Specify the extrusion direction 41 - -dm_plex_transform_extrude_thicknesses <t0,...,tl> - Specify thickness of each layer 42 43 Level: intermediate 44 45 .seealso: `DMExtrude()`, `DMPlexTransform`, `DMPlexTransformExtrudeSetThickness()`, `DMPlexTransformExtrudeSetTensor()` 46 @*/ 47 PetscErrorCode DMPlexExtrude(DM dm, PetscInt layers, PetscReal thickness, PetscBool tensor, PetscBool symmetric, const PetscReal normal[], const PetscReal thicknesses[], DM *edm) 48 { 49 DMPlexTransform tr; 50 DM cdm, ecdm; 51 const char *prefix; 52 PetscOptions options; 53 54 PetscFunctionBegin; 55 PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), &tr)); 56 PetscCall(DMPlexTransformSetDM(tr, dm)); 57 PetscCall(DMPlexTransformSetType(tr, DMPLEXEXTRUDE)); 58 PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix)); 59 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)tr, prefix)); 60 PetscCall(PetscObjectGetOptions((PetscObject)dm, &options)); 61 PetscCall(PetscObjectSetOptions((PetscObject)tr, options)); 62 PetscCall(DMPlexTransformExtrudeSetLayers(tr, layers)); 63 if (thickness > 0.) PetscCall(DMPlexTransformExtrudeSetThickness(tr, thickness)); 64 PetscCall(DMPlexTransformExtrudeSetTensor(tr, tensor)); 65 PetscCall(DMPlexTransformExtrudeSetSymmetric(tr, symmetric)); 66 if (normal) PetscCall(DMPlexTransformExtrudeSetNormal(tr, normal)); 67 if (thicknesses) PetscCall(DMPlexTransformExtrudeSetThicknesses(tr, layers, thicknesses)); 68 PetscCall(DMPlexTransformSetFromOptions(tr)); 69 PetscCall(PetscObjectSetOptions((PetscObject)tr, NULL)); 70 PetscCall(DMPlexTransformSetUp(tr)); 71 PetscCall(PetscObjectViewFromOptions((PetscObject)tr, NULL, "-dm_plex_transform_view")); 72 PetscCall(DMPlexTransformApply(tr, dm, edm)); 73 PetscCall(DMCopyDisc(dm, *edm)); 74 PetscCall(DMGetCoordinateDM(dm, &cdm)); 75 PetscCall(DMGetCoordinateDM(*edm, &ecdm)); 76 PetscCall(DMCopyDisc(cdm, ecdm)); 77 PetscCall(DMPlexTransformCreateDiscLabels(tr, *edm)); 78 PetscCall(DMPlexTransformDestroy(&tr)); 79 if (*edm) { 80 ((DM_Plex *)(*edm)->data)->printFEM = ((DM_Plex *)dm->data)->printFEM; 81 ((DM_Plex *)(*edm)->data)->printL2 = ((DM_Plex *)dm->data)->printL2; 82 } 83 PetscFunctionReturn(0); 84 } 85 86 PetscErrorCode DMExtrude_Plex(DM dm, PetscInt layers, DM *edm) 87 { 88 PetscFunctionBegin; 89 PetscCall(DMPlexExtrude(dm, layers, PETSC_DETERMINE, PETSC_TRUE, PETSC_FALSE, NULL, NULL, edm)); 90 PetscCall(DMViewFromOptions(*edm, NULL, "-check_extrude")); 91 PetscFunctionReturn(0); 92 } 93