1 #include <petsc/private/dmimpl.h> /*I "petscdm.h" I*/
2
3 PETSC_EXTERN PetscErrorCode DMIsForest(DM, PetscBool *);
4
5 DMGeneratorFunctionList DMGenerateList = NULL;
6 PetscBool DMGenerateRegisterAllCalled = PETSC_FALSE;
7
8 #if defined(PETSC_HAVE_TRIANGLE)
9 PETSC_EXTERN PetscErrorCode DMPlexGenerate_Triangle(DM, PetscBool, DM *);
10 PETSC_EXTERN PetscErrorCode DMPlexRefine_Triangle(DM, PetscReal *, DM *);
11 #endif
12 #if defined(PETSC_HAVE_TETGEN)
13 PETSC_EXTERN PetscErrorCode DMPlexGenerate_Tetgen(DM, PetscBool, DM *);
14 PETSC_EXTERN PetscErrorCode DMPlexRefine_Tetgen(DM, double *, DM *);
15 #endif
16 #if defined(PETSC_HAVE_CTETGEN)
17 PETSC_EXTERN PetscErrorCode DMPlexGenerate_CTetgen(DM, PetscBool, DM *);
18 PETSC_EXTERN PetscErrorCode DMPlexRefine_CTetgen(DM, double *, DM *);
19 #endif
20 #if defined(PETSC_HAVE_PRAGMATIC)
21 PETSC_EXTERN PetscErrorCode DMAdaptMetric_Pragmatic_Plex(DM, Vec, DMLabel, DMLabel, DM *);
22 #endif
23 #if defined(PETSC_HAVE_MMG)
24 PETSC_EXTERN PetscErrorCode DMAdaptMetric_Mmg_Plex(DM, Vec, DMLabel, DMLabel, DM *);
25 #endif
26 #if defined(PETSC_HAVE_PARMMG)
27 PETSC_EXTERN PetscErrorCode DMAdaptMetric_ParMmg_Plex(DM, Vec, DMLabel, DMLabel, DM *);
28 #endif
29 PETSC_EXTERN PetscErrorCode DMPlexTransformAdaptLabel(DM, Vec, DMLabel, DMLabel, DM *);
30 PETSC_EXTERN PetscErrorCode DMAdaptLabel_Forest(DM, Vec, DMLabel, DMLabel, DM *);
31
32 /*@C
33 DMGenerateRegisterAll - Registers all of the mesh generation methods in the `DM` package.
34
35 Not Collective
36
37 Level: advanced
38
39 .seealso: `DM`, `DMGenerateRegisterDestroy()`
40 @*/
DMGenerateRegisterAll(void)41 PetscErrorCode DMGenerateRegisterAll(void)
42 {
43 PetscFunctionBegin;
44 if (DMGenerateRegisterAllCalled) PetscFunctionReturn(PETSC_SUCCESS);
45 DMGenerateRegisterAllCalled = PETSC_TRUE;
46 #if defined(PETSC_HAVE_TRIANGLE)
47 PetscCall(DMGenerateRegister("triangle", DMPlexGenerate_Triangle, DMPlexRefine_Triangle, NULL, 1));
48 #endif
49 #if defined(PETSC_HAVE_CTETGEN)
50 PetscCall(DMGenerateRegister("ctetgen", DMPlexGenerate_CTetgen, DMPlexRefine_CTetgen, NULL, 2));
51 #endif
52 #if defined(PETSC_HAVE_TETGEN)
53 PetscCall(DMGenerateRegister("tetgen", DMPlexGenerate_Tetgen, DMPlexRefine_Tetgen, NULL, 2));
54 #endif
55 #if defined(PETSC_HAVE_PRAGMATIC)
56 PetscCall(DMGenerateRegister("pragmatic", NULL, NULL, DMAdaptMetric_Pragmatic_Plex, -1));
57 #endif
58 #if defined(PETSC_HAVE_MMG)
59 PetscCall(DMGenerateRegister("mmg", NULL, NULL, DMAdaptMetric_Mmg_Plex, -1));
60 #endif
61 #if defined(PETSC_HAVE_PARMMG)
62 PetscCall(DMGenerateRegister("parmmg", NULL, NULL, DMAdaptMetric_ParMmg_Plex, -1));
63 #endif
64 PetscCall(DMGenerateRegister("cellrefiner", NULL, NULL, DMPlexTransformAdaptLabel, -1));
65 PetscCall(DMGenerateRegister("forest", NULL, NULL, DMAdaptLabel_Forest, -1));
66 PetscFunctionReturn(PETSC_SUCCESS);
67 }
68
69 /*@C
70 DMGenerateRegister - Adds a grid generator to `DM`
71
72 Not Collective, No Fortran Support
73
74 Input Parameters:
75 + sname - name of a new user-defined grid generator
76 . fnc - generator function
77 . rfnc - refinement function
78 . alfnc - adapt by label function
79 - dim - dimension of boundary of domain
80
81 Example Usage:
82 .vb
83 DMGenerateRegister("my_generator", MyGeneratorCreate, MyGeneratorRefiner, MyGeneratorAdaptor, dim);
84 .ve
85
86 Then, your generator can be chosen with the procedural interface via
87 .vb
88 DMGenerate(dm, "my_generator",...)
89 .ve
90 or at runtime via the option
91 .vb
92 -dm_generator my_generator
93 .ve
94
95 Level: advanced
96
97 Note:
98 `DMGenerateRegister()` may be called multiple times to add several user-defined generators
99
100 .seealso: `DM`, `DMGenerateRegisterAll()`, `DMPlexGenerate()`, `DMGenerateRegisterDestroy()`
101 @*/
DMGenerateRegister(const char sname[],PetscErrorCode (* fnc)(DM,PetscBool,DM *),PetscErrorCode (* rfnc)(DM,PetscReal *,DM *),PetscErrorCode (* alfnc)(DM,Vec,DMLabel,DMLabel,DM *),PetscInt dim)102 PetscErrorCode DMGenerateRegister(const char sname[], PetscErrorCode (*fnc)(DM, PetscBool, DM *), PetscErrorCode (*rfnc)(DM, PetscReal *, DM *), PetscErrorCode (*alfnc)(DM, Vec, DMLabel, DMLabel, DM *), PetscInt dim)
103 {
104 DMGeneratorFunctionList entry;
105
106 PetscFunctionBegin;
107 PetscCall(PetscNew(&entry));
108 PetscCall(PetscStrallocpy(sname, &entry->name));
109 entry->generate = fnc;
110 entry->refine = rfnc;
111 entry->adapt = alfnc;
112 entry->dim = dim;
113 entry->next = NULL;
114 if (!DMGenerateList) DMGenerateList = entry;
115 else {
116 DMGeneratorFunctionList fl = DMGenerateList;
117 while (fl->next) fl = fl->next;
118 fl->next = entry;
119 }
120 PetscFunctionReturn(PETSC_SUCCESS);
121 }
122
DMGenerateRegisterDestroy(void)123 PetscErrorCode DMGenerateRegisterDestroy(void)
124 {
125 DMGeneratorFunctionList next, fl;
126
127 PetscFunctionBegin;
128 next = fl = DMGenerateList;
129 while (next) {
130 next = fl ? fl->next : NULL;
131 if (fl) PetscCall(PetscFree(fl->name));
132 PetscCall(PetscFree(fl));
133 fl = next;
134 }
135 DMGenerateList = NULL;
136 DMGenerateRegisterAllCalled = PETSC_FALSE;
137 PetscFunctionReturn(PETSC_SUCCESS);
138 }
139
140 /*@
141 DMAdaptLabel - Adapt a `DM` based on a `DMLabel` with values interpreted as coarsening and refining flags. Specific implementations of `DM` maybe have
142 specialized flags, but all implementations should accept flag values `DM_ADAPT_DETERMINE`, `DM_ADAPT_KEEP`, `DM_ADAPT_REFINE`, and,
143 `DM_ADAPT_COARSEN`.
144
145 Collective
146
147 Input Parameters:
148 + dm - the pre-adaptation `DM` object
149 - label - label with the flags
150
151 Output Parameters:
152 . dmAdapt - the adapted `DM` object: may be `NULL` if an adapted `DM` could not be produced.
153
154 Level: intermediate
155
156 .seealso: `DM`, `DMAdaptMetric()`, `DMCoarsen()`, `DMRefine()`
157 @*/
DMAdaptLabel(DM dm,DMLabel label,DM * dmAdapt)158 PetscErrorCode DMAdaptLabel(DM dm, DMLabel label, DM *dmAdapt)
159 {
160 DMGeneratorFunctionList fl;
161 char adaptname[PETSC_MAX_PATH_LEN];
162 const char *name;
163 PetscInt dim;
164 PetscBool flg, isForest, found = PETSC_FALSE;
165
166 PetscFunctionBegin;
167 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
168 if (label) PetscValidHeaderSpecific(label, DMLABEL_CLASSID, 2);
169 PetscAssertPointer(dmAdapt, 3);
170 *dmAdapt = NULL;
171 PetscCall(DMGetDimension(dm, &dim));
172 PetscCall(DMIsForest(dm, &isForest));
173 name = isForest ? "forest" : "cellrefiner";
174 PetscCall(PetscOptionsGetString(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dm_adaptor", adaptname, sizeof(adaptname), &flg));
175 if (flg) name = adaptname;
176
177 fl = DMGenerateList;
178 while (fl) {
179 PetscCall(PetscStrcmp(fl->name, name, &flg));
180 if (flg) {
181 PetscCall((*fl->adapt)(dm, NULL, label, NULL, dmAdapt));
182 found = PETSC_TRUE;
183 }
184 fl = fl->next;
185 }
186 PetscCheck(found, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Grid adaptor %s not registered; you may need to add --download-%s to your ./configure options", name, name);
187 if (*dmAdapt) {
188 PetscCtx ctx;
189
190 (*dmAdapt)->prealloc_only = dm->prealloc_only; /* maybe this should go .... */
191 PetscCall(PetscFree((*dmAdapt)->vectype));
192 PetscCall(PetscStrallocpy(dm->vectype, (char **)&(*dmAdapt)->vectype));
193 PetscCall(PetscFree((*dmAdapt)->mattype));
194 PetscCall(PetscStrallocpy(dm->mattype, (char **)&(*dmAdapt)->mattype));
195 PetscCall(DMGetApplicationContext(dm, &ctx));
196 PetscCall(DMSetApplicationContext(*dmAdapt, ctx));
197 }
198 PetscFunctionReturn(PETSC_SUCCESS);
199 }
200
201 /*@
202 DMAdaptMetric - Generates a mesh adapted to the specified metric field.
203
204 Input Parameters:
205 + dm - The DM object
206 . metric - The metric to which the mesh is adapted, defined vertex-wise.
207 . bdLabel - Label for boundary tags, which will be preserved in the output mesh. `bdLabel` should be `NULL` if there is no such label, and should be different from "_boundary_".
208 - rgLabel - Label for cell tags, which will be preserved in the output mesh. `rgLabel` should be `NULL` if there is no such label, and should be different from "_regions_".
209
210 Output Parameter:
211 . dmAdapt - Pointer to the `DM` object containing the adapted mesh
212
213 Note:
214 The label in the adapted mesh will be registered under the name of the input `DMLabel` object
215
216 Level: advanced
217
218 .seealso: `DMAdaptLabel()`, `DMCoarsen()`, `DMRefine()`
219 @*/
DMAdaptMetric(DM dm,Vec metric,DMLabel bdLabel,DMLabel rgLabel,DM * dmAdapt)220 PetscErrorCode DMAdaptMetric(DM dm, Vec metric, DMLabel bdLabel, DMLabel rgLabel, DM *dmAdapt)
221 {
222 DMGeneratorFunctionList fl;
223 char adaptname[PETSC_MAX_PATH_LEN];
224 const char *name;
225 const char *const adaptors[3] = {"pragmatic", "mmg", "parmmg"};
226 PetscInt dim;
227 PetscBool flg, found = PETSC_FALSE;
228
229 PetscFunctionBegin;
230 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
231 PetscValidHeaderSpecific(metric, VEC_CLASSID, 2);
232 if (bdLabel) PetscValidHeaderSpecific(bdLabel, DMLABEL_CLASSID, 3);
233 if (rgLabel) PetscValidHeaderSpecific(rgLabel, DMLABEL_CLASSID, 4);
234 PetscAssertPointer(dmAdapt, 5);
235 *dmAdapt = NULL;
236 PetscCall(DMGetDimension(dm, &dim));
237 PetscCall(PetscOptionsGetString(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dm_adaptor", adaptname, sizeof(adaptname), &flg));
238
239 /* Default to Mmg in serial and ParMmg in parallel */
240 if (flg) name = adaptname;
241 else {
242 MPI_Comm comm;
243 PetscMPIInt size;
244
245 PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
246 PetscCallMPI(MPI_Comm_size(comm, &size));
247 if (size == 1) name = adaptors[1];
248 else name = adaptors[2];
249 }
250
251 fl = DMGenerateList;
252 while (fl) {
253 PetscCall(PetscStrcmp(fl->name, name, &flg));
254 if (flg) {
255 PetscCall((*fl->adapt)(dm, metric, bdLabel, rgLabel, dmAdapt));
256 found = PETSC_TRUE;
257 }
258 fl = fl->next;
259 }
260 PetscCheck(found, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Grid adaptor %s not registered; you may need to add --download-%s to your ./configure options", name, name);
261 if (*dmAdapt) {
262 (*dmAdapt)->prealloc_only = dm->prealloc_only; /* maybe this should go .... */
263 PetscCall(PetscFree((*dmAdapt)->vectype));
264 PetscCall(PetscStrallocpy(dm->vectype, (char **)&(*dmAdapt)->vectype));
265 PetscCall(PetscFree((*dmAdapt)->mattype));
266 PetscCall(PetscStrallocpy(dm->mattype, (char **)&(*dmAdapt)->mattype));
267 }
268 PetscFunctionReturn(PETSC_SUCCESS);
269 }
270