1 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/
2 #include <petscdmplextransform.h>
3
4 /*@
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 product cells
12 . symmetric - Flag to extrude symmetrically about the surface
13 . periodic - Flag to extrude periodically
14 . normal - Surface normal vector, or `NULL`
15 . thicknesses - Thickness of each layer, or `NULL`
16 - activeLabel - `DMLabel` to extrude from, or `NULL` to extrude entire mesh
17
18 Output Parameter:
19 . edm - The volumetric mesh
20
21 Options Database Keys:
22 + -dm_plex_transform_extrude_thickness <t> - The total thickness of extruded layers
23 . -dm_plex_transform_extrude_use_tensor <bool> - Use tensor cells when extruding
24 . -dm_plex_transform_extrude_symmetric <bool> - Extrude layers symmetrically about the surface
25 . -dm_plex_transform_extrude_periodic <bool> - Extrude layers periodically
26 . -dm_plex_transform_extrude_normal <n0,...,nd> - Specify the extrusion direction
27 - -dm_plex_transform_extrude_thicknesses <t0,...,tl> - Specify thickness of each layer
28
29 Level: intermediate
30
31 Notes:
32 Extrusion is implemented as a `DMPlexTransform`, so that new mesh points are produced from old mesh points. In the example below,
33 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,
34 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
35 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
36 would v1, v2, e1 and e2.
37
38 .vb
39 v2----- e6 -----v5
40 | |
41 e2 face2 e4
42 | |
43 v1----- e5 -----v4
44 | |
45 e1 face1 e3
46 | |
47 v0--- original ----v3
48 .ve
49
50 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMExtrude()`, `DMPlexTransform`, `DMPlexTransformExtrudeSetThickness()`, `DMPlexTransformExtrudeSetTensor()`
51 @*/
DMPlexExtrude(DM dm,PetscInt layers,PetscReal thickness,PetscBool tensor,PetscBool symmetric,PetscBool periodic,const PetscReal normal[],const PetscReal thicknesses[],DMLabel activeLabel,DM * edm)52 PetscErrorCode DMPlexExtrude(DM dm, PetscInt layers, PetscReal thickness, PetscBool tensor, PetscBool symmetric, PetscBool periodic, const PetscReal normal[], const PetscReal thicknesses[], DMLabel activeLabel, DM *edm)
53 {
54 DMPlexTransform tr;
55 DM cdm;
56 PetscObject disc;
57 PetscClassId id;
58 const char *prefix;
59 PetscOptions options;
60 PetscBool cutMarker = PETSC_FALSE;
61
62 PetscFunctionBegin;
63 PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), &tr));
64 PetscCall(PetscObjectSetName((PetscObject)tr, "Extrusion Transform"));
65 PetscCall(DMPlexTransformSetDM(tr, dm));
66 PetscCall(DMPlexTransformSetType(tr, DMPLEXEXTRUDETYPE));
67 PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix));
68 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)tr, prefix));
69 PetscCall(PetscObjectGetOptions((PetscObject)dm, &options));
70 PetscCall(PetscObjectSetOptions((PetscObject)tr, options));
71 if (activeLabel) PetscCall(DMPlexTransformSetActive(tr, activeLabel));
72 PetscCall(DMPlexTransformExtrudeSetLayers(tr, layers));
73 if (thickness > 0.) PetscCall(DMPlexTransformExtrudeSetThickness(tr, thickness));
74 PetscCall(DMPlexTransformExtrudeSetTensor(tr, tensor));
75 PetscCall(DMPlexTransformExtrudeSetSymmetric(tr, symmetric));
76 PetscCall(DMPlexTransformExtrudeSetPeriodic(tr, periodic));
77 if (normal) PetscCall(DMPlexTransformExtrudeSetNormal(tr, normal));
78 if (thicknesses) PetscCall(DMPlexTransformExtrudeSetThicknesses(tr, layers, thicknesses));
79 PetscCall(DMPlexTransformSetFromOptions(tr));
80 PetscCall(PetscObjectSetOptions((PetscObject)tr, NULL));
81 PetscCall(DMPlexTransformSetUp(tr));
82 PetscCall(PetscObjectViewFromOptions((PetscObject)tr, NULL, "-dm_plex_transform_view"));
83 PetscCall(DMPlexTransformApply(tr, dm, edm));
84 PetscCall(DMCopyDisc(dm, *edm));
85 // Handle periodic viewing
86 PetscCall(PetscOptionsGetBool(options, ((PetscObject)dm)->prefix, "-dm_plex_periodic_cut", &cutMarker, NULL));
87 PetscCall(DMPlexTransformExtrudeGetPeriodic(tr, &periodic));
88 if (periodic && cutMarker) {
89 DMLabel cutLabel;
90 PetscInt dim, pStart, pEnd;
91
92 PetscCall(DMGetDimension(dm, &dim));
93 PetscCall(DMCreateLabel(*edm, "periodic_cut"));
94 PetscCall(DMGetLabel(*edm, "periodic_cut", &cutLabel));
95 PetscCall(DMPlexGetChart(dm, &pStart, &pEnd));
96 for (PetscInt p = pStart; p < pEnd; ++p) {
97 DMPolytopeType ct;
98 DMPolytopeType *rct;
99 PetscInt *rsize, *rcone, *rornt;
100 PetscInt Nct;
101
102 PetscCall(DMPlexGetCellType(dm, p, &ct));
103 PetscCall(DMPlexTransformCellTransform(tr, ct, p, NULL, &Nct, &rct, &rsize, &rcone, &rornt));
104 for (PetscInt n = 0; n < Nct; ++n) {
105 PetscInt pNew, pdim = DMPolytopeTypeGetDim(rct[n]);
106
107 if (ct == rct[n] || pdim > dim) {
108 PetscCall(DMPlexTransformGetTargetPoint(tr, ct, rct[n], p, 0, &pNew));
109 PetscCall(DMLabelSetValue(cutLabel, pNew, !pdim ? 1 : 2));
110 }
111 }
112 }
113 }
114 // It is too hard to raise the dimension of a discretization, so just remake it
115 PetscCall(DMGetCoordinateDM(dm, &cdm));
116 PetscCall(DMGetField(cdm, 0, NULL, &disc));
117 PetscCall(PetscObjectGetClassId(disc, &id));
118 if (id == PETSCFE_CLASSID) {
119 PetscSpace sp;
120 PetscInt deg;
121
122 PetscCall(PetscFEGetBasisSpace((PetscFE)disc, &sp));
123 PetscCall(PetscSpaceGetDegree(sp, °, NULL));
124 PetscCall(DMPlexCreateCoordinateSpace(*edm, deg, PETSC_FALSE, PETSC_TRUE));
125 }
126 PetscCall(DMPlexTransformCreateDiscLabels(tr, *edm));
127 PetscCall(DMPlexTransformDestroy(&tr));
128 PetscCall(DMPlexCopy_Internal(dm, PETSC_FALSE, PETSC_FALSE, *edm));
129 PetscFunctionReturn(PETSC_SUCCESS);
130 }
131
DMExtrude_Plex(DM dm,PetscInt layers,DM * edm)132 PetscErrorCode DMExtrude_Plex(DM dm, PetscInt layers, DM *edm)
133 {
134 PetscFunctionBegin;
135 PetscCall(DMPlexExtrude(dm, layers, PETSC_DETERMINE, PETSC_TRUE, PETSC_FALSE, PETSC_FALSE, NULL, NULL, NULL, edm));
136 PetscCall(DMSetMatType(*edm, dm->mattype));
137 PetscCall(DMViewFromOptions(*edm, NULL, "-check_extrude"));
138 PetscFunctionReturn(PETSC_SUCCESS);
139 }
140