Actual source code: plexextrude.c
1: #include <petsc/private/dmpleximpl.h>
2: #include <petscdmplextransform.h>
4: /*@C
5: DMPlexExtrude - Extrude a volumetric mesh from the input surface mesh
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
16: Output Parameter:
17: . edm - The volumetric mesh
19: Options Database Keys:
20: + -dm_plex_transform_extrude_thickness <t> - The total thickness of extruded layers
21: . -dm_plex_transform_extrude_use_tensor <bool> - Use tensor cells when extruding
22: . -dm_plex_transform_extrude_symmetric <bool> - Extrude layers symmetrically about the surface
23: . -dm_plex_transform_extrude_normal <n0,...,nd> - Specify the extrusion direction
24: - -dm_plex_transform_extrude_thicknesses <t0,...,tl> - Specify thickness of each layer
26: Level: intermediate
28: Notes:
29: Extrusion is implemented as a `DMPlexTransform`, so that new mesh points are produced from old mesh points. In the example below,
30: 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,
31: 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
32: 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
33: would v1, v2, e1 and e2.
35: .vb
36: v2----- e6 -----v5
37: | |
38: e2 face2 e4
39: | |
40: v1----- e5 -----v4
41: | |
42: e1 face1 e3
43: | |
44: v0--- original ----v3
45: .ve
47: .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMExtrude()`, `DMPlexTransform`, `DMPlexTransformExtrudeSetThickness()`, `DMPlexTransformExtrudeSetTensor()`
48: @*/
49: PetscErrorCode DMPlexExtrude(DM dm, PetscInt layers, PetscReal thickness, PetscBool tensor, PetscBool symmetric, const PetscReal normal[], const PetscReal thicknesses[], DM *edm)
50: {
51: DMPlexTransform tr;
52: DM cdm;
53: PetscObject disc;
54: PetscClassId id;
55: const char *prefix;
56: PetscOptions options;
58: PetscFunctionBegin;
59: PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), &tr));
60: PetscCall(DMPlexTransformSetDM(tr, dm));
61: PetscCall(DMPlexTransformSetType(tr, DMPLEXEXTRUDE));
62: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix));
63: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)tr, prefix));
64: PetscCall(PetscObjectGetOptions((PetscObject)dm, &options));
65: PetscCall(PetscObjectSetOptions((PetscObject)tr, options));
66: PetscCall(DMPlexTransformExtrudeSetLayers(tr, layers));
67: if (thickness > 0.) PetscCall(DMPlexTransformExtrudeSetThickness(tr, thickness));
68: PetscCall(DMPlexTransformExtrudeSetTensor(tr, tensor));
69: PetscCall(DMPlexTransformExtrudeSetSymmetric(tr, symmetric));
70: if (normal) PetscCall(DMPlexTransformExtrudeSetNormal(tr, normal));
71: if (thicknesses) PetscCall(DMPlexTransformExtrudeSetThicknesses(tr, layers, thicknesses));
72: PetscCall(DMPlexTransformSetFromOptions(tr));
73: PetscCall(PetscObjectSetOptions((PetscObject)tr, NULL));
74: PetscCall(DMPlexTransformSetUp(tr));
75: PetscCall(PetscObjectViewFromOptions((PetscObject)tr, NULL, "-dm_plex_transform_view"));
76: PetscCall(DMPlexTransformApply(tr, dm, edm));
77: PetscCall(DMCopyDisc(dm, *edm));
78: // It is too hard to raise the dimension of a discretization, so just remake it
79: PetscCall(DMGetCoordinateDM(dm, &cdm));
80: PetscCall(DMGetField(cdm, 0, NULL, &disc));
81: PetscCall(PetscObjectGetClassId(disc, &id));
82: if (id == PETSCFE_CLASSID) {
83: PetscSpace sp;
84: PetscInt deg;
86: PetscCall(PetscFEGetBasisSpace((PetscFE)disc, &sp));
87: PetscCall(PetscSpaceGetDegree(sp, °, NULL));
88: PetscCall(DMPlexCreateCoordinateSpace(*edm, deg, NULL));
89: }
90: PetscCall(DMPlexTransformCreateDiscLabels(tr, *edm));
91: PetscCall(DMPlexTransformDestroy(&tr));
92: if (*edm) {
93: ((DM_Plex *)(*edm)->data)->printFEM = ((DM_Plex *)dm->data)->printFEM;
94: ((DM_Plex *)(*edm)->data)->printL2 = ((DM_Plex *)dm->data)->printL2;
95: }
96: PetscFunctionReturn(PETSC_SUCCESS);
97: }
99: PetscErrorCode DMExtrude_Plex(DM dm, PetscInt layers, DM *edm)
100: {
101: PetscFunctionBegin;
102: PetscCall(DMPlexExtrude(dm, layers, PETSC_DETERMINE, PETSC_TRUE, PETSC_FALSE, NULL, NULL, edm));
103: PetscCall(DMViewFromOptions(*edm, NULL, "-check_extrude"));
104: PetscFunctionReturn(PETSC_SUCCESS);
105: }