Actual source code: plexgenerate.c
1: #include <petsc/private/dmpleximpl.h>
3: /*@C
4: DMPlexInvertCell - Flips cell orientations since `DMPLEX` stores some of them internally with outward normals.
6: Input Parameters:
7: + cellType - The cell type
8: - cone - The incoming cone
10: Output Parameter:
11: . cone - The inverted cone (in-place)
13: Level: developer
15: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMPolytopeType`, `DMPlexGenerate()`
16: @*/
17: PetscErrorCode DMPlexInvertCell(DMPolytopeType cellType, PetscInt cone[])
18: {
19: #define SWAPCONE(cone, i, j) \
20: do { \
21: PetscInt _cone_tmp; \
22: _cone_tmp = (cone)[i]; \
23: (cone)[i] = (cone)[j]; \
24: (cone)[j] = _cone_tmp; \
25: } while (0)
27: switch (cellType) {
28: case DM_POLYTOPE_POINT:
29: break;
30: case DM_POLYTOPE_SEGMENT:
31: break;
32: case DM_POLYTOPE_POINT_PRISM_TENSOR:
33: break;
34: case DM_POLYTOPE_TRIANGLE:
35: break;
36: case DM_POLYTOPE_QUADRILATERAL:
37: break;
38: case DM_POLYTOPE_SEG_PRISM_TENSOR:
39: SWAPCONE(cone, 2, 3);
40: break;
41: case DM_POLYTOPE_TETRAHEDRON:
42: SWAPCONE(cone, 0, 1);
43: break;
44: case DM_POLYTOPE_HEXAHEDRON:
45: SWAPCONE(cone, 1, 3);
46: break;
47: case DM_POLYTOPE_TRI_PRISM:
48: SWAPCONE(cone, 1, 2);
49: break;
50: case DM_POLYTOPE_TRI_PRISM_TENSOR:
51: break;
52: case DM_POLYTOPE_QUAD_PRISM_TENSOR:
53: break;
54: case DM_POLYTOPE_PYRAMID:
55: SWAPCONE(cone, 1, 3);
56: break;
57: default:
58: break;
59: }
60: return 0;
61: #undef SWAPCONE
62: }
64: /*@C
65: DMPlexReorderCell - Flips cell orientations since `DMPLEX` stores some of them internally with outward normals.
67: Input Parameters:
68: + dm - The `DMPLEX` object
69: . cell - The cell
70: - cone - The incoming cone
72: Output Parameter:
73: . cone - The reordered cone (in-place)
75: Level: developer
77: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMPolytopeType`, `DMPlexGenerate()`
78: @*/
79: PetscErrorCode DMPlexReorderCell(DM dm, PetscInt cell, PetscInt cone[])
80: {
81: DMPolytopeType cellType;
83: DMPlexGetCellType(dm, cell, &cellType);
84: DMPlexInvertCell(cellType, cone);
85: return 0;
86: }
88: /*@C
89: DMPlexTriangleSetOptions - Set the options used for the Triangle mesh generator
91: Not Collective
93: Inputs Parameters:
94: + dm - The `DMPLEX` object
95: - opts - The command line options
97: Level: developer
99: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMPlexTetgenSetOptions()`, `DMPlexGenerate()`
100: @*/
101: PetscErrorCode DMPlexTriangleSetOptions(DM dm, const char *opts)
102: {
103: DM_Plex *mesh = (DM_Plex *)dm->data;
107: PetscFree(mesh->triangleOpts);
108: PetscStrallocpy(opts, &mesh->triangleOpts);
109: return 0;
110: }
112: /*@C
113: DMPlexTetgenSetOptions - Set the options used for the Tetgen mesh generator
115: Not Collective
117: Inputs Parameters:
118: + dm - The `DMPLEX` object
119: - opts - The command line options
121: Level: developer
123: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMPlexTriangleSetOptions()`, `DMPlexGenerate()`
124: @*/
125: PetscErrorCode DMPlexTetgenSetOptions(DM dm, const char *opts)
126: {
127: DM_Plex *mesh = (DM_Plex *)dm->data;
131: PetscFree(mesh->tetgenOpts);
132: PetscStrallocpy(opts, &mesh->tetgenOpts);
133: return 0;
134: }
136: /*@C
137: DMPlexGenerate - Generates a mesh.
139: Not Collective
141: Input Parameters:
142: + boundary - The `DMPLEX` boundary object
143: . name - The mesh generation package name
144: - interpolate - Flag to create intermediate mesh elements
146: Output Parameter:
147: . mesh - The `DMPLEX` object
149: Options Database Keys:
150: + -dm_plex_generate <name> - package to generate mesh, for example, triangle, ctetgen or tetgen
151: - -dm_generator <name> - package to generate mesh, for example, triangle, ctetgen or tetgen
153: Level: intermediate
155: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMPlexCreate()`, `DMRefine()`
156: @*/
157: PetscErrorCode DMPlexGenerate(DM boundary, const char name[], PetscBool interpolate, DM *mesh)
158: {
159: DMGeneratorFunctionList fl;
160: char genname[PETSC_MAX_PATH_LEN];
161: const char *suggestions;
162: PetscInt dim;
163: PetscBool flg;
167: DMGetDimension(boundary, &dim);
168: PetscOptionsGetString(((PetscObject)boundary)->options, ((PetscObject)boundary)->prefix, "-dm_generator", genname, sizeof(genname), &flg);
169: if (flg) name = genname;
170: else {
171: PetscOptionsGetString(((PetscObject)boundary)->options, ((PetscObject)boundary)->prefix, "-dm_plex_generate", genname, sizeof(genname), &flg);
172: if (flg) name = genname;
173: }
175: fl = DMGenerateList;
176: if (name) {
177: while (fl) {
178: PetscStrcmp(fl->name, name, &flg);
179: if (flg) {
180: (*fl->generate)(boundary, interpolate, mesh);
181: return 0;
182: }
183: fl = fl->next;
184: }
185: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Grid generator %s not registered; you may need to add --download-%s to your ./configure options", name, name);
186: } else {
187: while (fl) {
188: if (boundary->dim == fl->dim) {
189: (*fl->generate)(boundary, interpolate, mesh);
190: return 0;
191: }
192: fl = fl->next;
193: }
194: suggestions = "";
195: if (boundary->dim + 1 == 2) suggestions = " You may need to add --download-triangle to your ./configure options";
196: else if (boundary->dim + 1 == 3) suggestions = " You may need to add --download-ctetgen or --download-tetgen in your ./configure options";
197: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No grid generator of dimension %" PetscInt_FMT " registered%s", boundary->dim + 1, suggestions);
198: }
199: }