Actual source code: plexgenerate.c
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: }
72: /*@C
73: DMPlexTriangleSetOptions - Set the options used for the Triangle mesh generator
75: Not Collective
77: Inputs Parameters:
78: + dm - The DMPlex object
79: - opts - The command line options
81: Level: developer
83: .seealso: DMPlexTetgenSetOptions(), DMPlexGenerate()
84: @*/
85: PetscErrorCode DMPlexTriangleSetOptions(DM dm, const char *opts)
86: {
87: DM_Plex *mesh = (DM_Plex*) dm->data;
93: PetscFree(mesh->triangleOpts);
94: PetscStrallocpy(opts, &mesh->triangleOpts);
95: return(0);
96: }
98: /*@C
99: DMPlexTetgenSetOptions - Set the options used for the Tetgen mesh generator
101: Not Collective
103: Inputs Parameters:
104: + dm - The DMPlex object
105: - opts - The command line options
107: Level: developer
109: .seealso: DMPlexTriangleSetOptions(), DMPlexGenerate()
110: @*/
111: PetscErrorCode DMPlexTetgenSetOptions(DM dm, const char *opts)
112: {
113: DM_Plex *mesh = (DM_Plex*) dm->data;
119: PetscFree(mesh->tetgenOpts);
120: PetscStrallocpy(opts, &mesh->tetgenOpts);
121: return(0);
122: }
124: /*
125: Contains the list of registered DMPlexGenerators routines
126: */
127: PlexGeneratorFunctionList DMPlexGenerateList = NULL;
129: /*@C
130: DMPlexGenerate - Generates a mesh.
132: Not Collective
134: Input Parameters:
135: + boundary - The DMPlex boundary object
136: . name - The mesh generation package name
137: - interpolate - Flag to create intermediate mesh elements
139: Output Parameter:
140: . mesh - The DMPlex object
142: Options Database:
143: + -dm_plex_generate <name> - package to generate mesh, for example, triangle, ctetgen or tetgen
144: - -dm_plex_generator <name> - package to generate mesh, for example, triangle, ctetgen or tetgen (deprecated)
146: Level: intermediate
148: .seealso: DMPlexCreate(), DMRefine()
149: @*/
150: PetscErrorCode DMPlexGenerate(DM boundary, const char name[], PetscBool interpolate, DM *mesh)
151: {
152: PlexGeneratorFunctionList fl;
153: char genname[PETSC_MAX_PATH_LEN];
154: const char *suggestions;
155: PetscInt dim;
156: PetscBool flg;
157: PetscErrorCode ierr;
162: DMGetDimension(boundary, &dim);
163: PetscOptionsGetString(((PetscObject) boundary)->options,((PetscObject) boundary)->prefix, "-dm_plex_generator", genname, sizeof(genname), &flg);
164: if (flg) name = genname;
165: else {
166: PetscOptionsGetString(((PetscObject) boundary)->options,((PetscObject) boundary)->prefix, "-dm_plex_generate", genname, sizeof(genname), &flg);
167: if (flg) name = genname;
168: }
170: fl = DMPlexGenerateList;
171: if (name) {
172: while (fl) {
173: PetscStrcmp(fl->name,name,&flg);
174: if (flg) {
175: (*fl->generate)(boundary,interpolate,mesh);
176: return(0);
177: }
178: fl = fl->next;
179: }
180: 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);
181: } else {
182: while (fl) {
183: if (boundary->dim == fl->dim) {
184: (*fl->generate)(boundary,interpolate,mesh);
185: return(0);
186: }
187: fl = fl->next;
188: }
189: suggestions = "";
190: if (boundary->dim+1 == 2) suggestions = " You may need to add --download-triangle to your ./configure options";
191: else if (boundary->dim+1 == 3) suggestions = " You may need to add --download-ctetgen or --download-tetgen in your ./configure options";
192: SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"No grid generator of dimension %D registered%s",boundary->dim+1,suggestions);
193: }
194: }
196: /*@C
197: DMPlexGenerateRegister - Adds a grid generator to DMPlex
199: Not Collective
201: Input Parameters:
202: + name_solver - name of a new user-defined grid generator
203: . fnc - generator function
204: . rfnc - refinement function
205: . alfnc - adapt by label function
206: - dim - dimension of boundary of domain
208: Notes:
209: DMPlexGenerateRegister() may be called multiple times to add several user-defined solvers.
211: Sample usage:
212: .vb
213: DMPlexGenerateRegister("my_generator",MyGeneratorCreate,MyGeneratorRefiner,MyGeneratorAdaptor,dim);
214: .ve
216: Then, your generator can be chosen with the procedural interface via
217: $ DMPlexGenerate(dm,"my_generator",...)
218: or at runtime via the option
219: $ -dm_plex_generator my_generator
221: Level: advanced
223: .seealso: DMPlexGenerateRegisterAll(), DMPlexGenerate(), DMPlexGenerateRegisterDestroy()
225: @*/
226: PetscErrorCode DMPlexGenerateRegister(const char sname[], PetscErrorCode (*fnc)(DM, PetscBool, DM*), PetscErrorCode (*rfnc)(DM, PetscReal*, DM*), PetscErrorCode (*alfnc)(DM, DMLabel, DM*), PetscInt dim)
227: {
228: PlexGeneratorFunctionList entry;
229: PetscErrorCode ierr;
232: PetscNew(&entry);
233: PetscStrallocpy(sname,&entry->name);
234: entry->generate = fnc;
235: entry->refine = rfnc;
236: entry->adaptlabel = alfnc;
237: entry->dim = dim;
238: entry->next = NULL;
239: if (!DMPlexGenerateList) DMPlexGenerateList = entry;
240: else {
241: PlexGeneratorFunctionList fl = DMPlexGenerateList;
242: while (fl->next) fl = fl->next;
243: fl->next = entry;
244: }
245: return(0);
246: }
248: extern PetscBool DMPlexGenerateRegisterAllCalled;
250: PetscErrorCode DMPlexGenerateRegisterDestroy(void)
251: {
252: PlexGeneratorFunctionList next,fl;
253: PetscErrorCode ierr;
256: next = fl = DMPlexGenerateList;
257: while (next) {
258: next = fl ? fl->next : NULL;
259: if (fl) {PetscFree(fl->name);}
260: PetscFree(fl);
261: fl = next;
262: }
263: DMPlexGenerateList = NULL;
264: DMPlexGenerateRegisterAllCalled = PETSC_FALSE;
265: return(0);
266: }