Actual source code: plexgenerate.c
petsc-master 2021-02-25
1: #include <petsc/private/dmpleximpl.h>
3: /*@C
4: DMPlexInvertCell - Flips cell orientations since Plex 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: 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)
28: switch (cellType) {
29: case DM_POLYTOPE_POINT: break;
30: case DM_POLYTOPE_SEGMENT: break;
31: case DM_POLYTOPE_POINT_PRISM_TENSOR: break;
32: case DM_POLYTOPE_TRIANGLE: break;
33: case DM_POLYTOPE_QUADRILATERAL: break;
34: case DM_POLYTOPE_SEG_PRISM_TENSOR: SWAPCONE(cone,2,3); break;
35: case DM_POLYTOPE_TETRAHEDRON: SWAPCONE(cone,0,1); break;
36: case DM_POLYTOPE_HEXAHEDRON: SWAPCONE(cone,1,3); break;
37: case DM_POLYTOPE_TRI_PRISM: SWAPCONE(cone,1,2); break;
38: case DM_POLYTOPE_TRI_PRISM_TENSOR: break;
39: case DM_POLYTOPE_QUAD_PRISM_TENSOR: break;
40: default: break;
41: }
42: return(0);
43: #undef SWAPCONE
44: }
46: /*@C
47: DMPlexReorderCell - Flips cell orientations since Plex stores some of them internally with outward normals.
49: Input Parameters:
50: + dm - The DMPlex object
51: . cell - The cell
52: - cone - The incoming cone
54: Output Parameter:
55: . cone - The reordered cone (in-place)
57: Level: developer
59: .seealso: DMPlexGenerate()
60: @*/
61: PetscErrorCode DMPlexReorderCell(DM dm, PetscInt cell, PetscInt cone[])
62: {
63: DMPolytopeType cellType;
67: DMPlexGetCellType(dm, cell, &cellType);
68: DMPlexInvertCell(cellType, cone);
69: return(0);
70: }
73: /*@C
74: DMPlexTriangleSetOptions - Set the options used for the Triangle mesh generator
76: Not Collective
78: Inputs Parameters:
79: + dm - The DMPlex object
80: - opts - The command line options
82: Level: developer
84: .seealso: DMPlexTetgenSetOptions(), DMPlexGenerate()
85: @*/
86: PetscErrorCode DMPlexTriangleSetOptions(DM dm, const char *opts)
87: {
88: DM_Plex *mesh = (DM_Plex*) dm->data;
94: PetscFree(mesh->triangleOpts);
95: PetscStrallocpy(opts, &mesh->triangleOpts);
96: return(0);
97: }
99: /*@C
100: DMPlexTetgenSetOptions - Set the options used for the Tetgen mesh generator
102: Not Collective
104: Inputs Parameters:
105: + dm - The DMPlex object
106: - opts - The command line options
108: Level: developer
110: .seealso: DMPlexTriangleSetOptions(), DMPlexGenerate()
111: @*/
112: PetscErrorCode DMPlexTetgenSetOptions(DM dm, const char *opts)
113: {
114: DM_Plex *mesh = (DM_Plex*) dm->data;
120: PetscFree(mesh->tetgenOpts);
121: PetscStrallocpy(opts, &mesh->tetgenOpts);
122: return(0);
123: }
125: /*
126: Contains the list of registered DMPlexGenerators routines
127: */
128: PlexGeneratorFunctionList DMPlexGenerateList = NULL;
130: /*@C
131: DMPlexGenerate - Generates a mesh.
133: Not Collective
135: Input Parameters:
136: + boundary - The DMPlex boundary object
137: . name - The mesh generation package name
138: - interpolate - Flag to create intermediate mesh elements
140: Output Parameter:
141: . mesh - The DMPlex object
143: Options Database:
144: + -dm_plex_generate <name> - package to generate mesh, for example, triangle, ctetgen or tetgen
145: - -dm_plex_generator <name> - package to generate mesh, for example, triangle, ctetgen or tetgen (deprecated)
147: Level: intermediate
149: .seealso: DMPlexCreate(), DMRefine()
150: @*/
151: PetscErrorCode DMPlexGenerate(DM boundary, const char name[], PetscBool interpolate, DM *mesh)
152: {
153: PlexGeneratorFunctionList fl;
154: char genname[PETSC_MAX_PATH_LEN];
155: const char *suggestions;
156: PetscInt dim;
157: PetscBool flg;
158: PetscErrorCode ierr;
163: DMGetDimension(boundary, &dim);
164: PetscOptionsGetString(((PetscObject) boundary)->options,((PetscObject) boundary)->prefix, "-dm_plex_generator", genname, sizeof(genname), &flg);
165: if (flg) name = genname;
166: else {
167: PetscOptionsGetString(((PetscObject) boundary)->options,((PetscObject) boundary)->prefix, "-dm_plex_generate", genname, sizeof(genname), &flg);
168: if (flg) name = genname;
169: }
171: fl = DMPlexGenerateList;
172: if (name) {
173: while (fl) {
174: PetscStrcmp(fl->name,name,&flg);
175: if (flg) {
176: (*fl->generate)(boundary,interpolate,mesh);
177: return(0);
178: }
179: fl = fl->next;
180: }
181: SETERRQ2(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);
182: } else {
183: while (fl) {
184: if (boundary->dim == fl->dim) {
185: (*fl->generate)(boundary,interpolate,mesh);
186: return(0);
187: }
188: fl = fl->next;
189: }
190: suggestions = "";
191: if (boundary->dim+1 == 2) suggestions = " You may need to add --download-triangle to your ./configure options";
192: else if (boundary->dim+1 == 3) suggestions = " You may need to add --download-ctetgen or --download-tetgen in your ./configure options";
193: SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"No grid generator of dimension %D registered%s",boundary->dim+1,suggestions);
194: }
195: }
197: /*@C
198: DMPlexGenerateRegister - Adds a grid generator to DMPlex
200: Not Collective
202: Input Parameters:
203: + name_solver - name of a new user-defined grid generator
204: . fnc - generator function
205: . rfnc - refinement function
206: . alfnc - adapt by label function
207: - dim - dimension of boundary of domain
209: Notes:
210: DMPlexGenerateRegister() may be called multiple times to add several user-defined solvers.
212: Sample usage:
213: .vb
214: DMPlexGenerateRegister("my_generator",MyGeneratorCreate,MyGeneratorRefiner,MyGeneratorAdaptor,dim);
215: .ve
217: Then, your generator can be chosen with the procedural interface via
218: $ DMPlexGenerate(dm,"my_generator",...)
219: or at runtime via the option
220: $ -dm_plex_generator my_generator
222: Level: advanced
224: .seealso: DMPlexGenerateRegisterAll(), DMPlexGenerate(), DMPlexGenerateRegisterDestroy()
226: @*/
227: PetscErrorCode DMPlexGenerateRegister(const char sname[], PetscErrorCode (*fnc)(DM, PetscBool, DM*), PetscErrorCode (*rfnc)(DM, PetscReal*, DM*), PetscErrorCode (*alfnc)(DM, DMLabel, DM*), PetscInt dim)
228: {
229: PlexGeneratorFunctionList entry;
230: PetscErrorCode ierr;
233: PetscNew(&entry);
234: PetscStrallocpy(sname,&entry->name);
235: entry->generate = fnc;
236: entry->refine = rfnc;
237: entry->adaptlabel = alfnc;
238: entry->dim = dim;
239: entry->next = NULL;
240: if (!DMPlexGenerateList) DMPlexGenerateList = entry;
241: else {
242: PlexGeneratorFunctionList fl = DMPlexGenerateList;
243: while (fl->next) fl = fl->next;
244: fl->next = entry;
245: }
246: return(0);
247: }
249: extern PetscBool DMPlexGenerateRegisterAllCalled;
251: PetscErrorCode DMPlexGenerateRegisterDestroy(void)
252: {
253: PlexGeneratorFunctionList next,fl;
254: PetscErrorCode ierr;
257: next = fl = DMPlexGenerateList;
258: while (next) {
259: next = fl ? fl->next : NULL;
260: if (fl) {PetscFree(fl->name);}
261: PetscFree(fl);
262: fl = next;
263: }
264: DMPlexGenerateList = NULL;
265: DMPlexGenerateRegisterAllCalled = PETSC_FALSE;
266: return(0);
267: }