1 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/
2 #include <petsc/private/petscfeimpl.h> /* For PetscFEInterpolate_Static() */
3
4 #include <petscdmplextransform.h> /*I "petscdmplextransform.h" I*/
5 #include <petscsf.h>
6
7 /*@
8 DMPlexCreateProcessSF - Create an `PetscSF` which just has process connectivity
9
10 Collective
11
12 Input Parameters:
13 + dm - The `DM`
14 - sfPoint - The `PetscSF` which encodes point connectivity
15
16 Output Parameters:
17 + processRanks - A list of process neighbors, or `NULL`
18 - sfProcess - An `PetscSF` encoding the process connectivity, or `NULL`
19
20 Level: developer
21
22 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `PetscSF`, `PetscSFCreate()`, `DMPlexCreateTwoSidedProcessSF()`
23 @*/
DMPlexCreateProcessSF(DM dm,PetscSF sfPoint,IS * processRanks,PetscSF * sfProcess)24 PetscErrorCode DMPlexCreateProcessSF(DM dm, PetscSF sfPoint, IS *processRanks, PetscSF *sfProcess)
25 {
26 PetscInt numRoots, numLeaves, l;
27 const PetscInt *localPoints;
28 const PetscSFNode *remotePoints;
29 PetscInt *localPointsNew;
30 PetscSFNode *remotePointsNew;
31 PetscMPIInt *ranks;
32 PetscInt *ranksNew;
33 PetscMPIInt size;
34
35 PetscFunctionBegin;
36 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
37 PetscValidHeaderSpecific(sfPoint, PETSCSF_CLASSID, 2);
38 if (processRanks) PetscAssertPointer(processRanks, 3);
39 if (sfProcess) PetscAssertPointer(sfProcess, 4);
40 PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)dm), &size));
41 PetscCall(PetscSFGetGraph(sfPoint, &numRoots, &numLeaves, &localPoints, &remotePoints));
42 PetscCall(PetscMalloc1(numLeaves, &ranks));
43 for (l = 0; l < numLeaves; ++l) ranks[l] = (PetscMPIInt)remotePoints[l].rank;
44 PetscCall(PetscSortRemoveDupsMPIInt(&numLeaves, ranks));
45 PetscCall(PetscMalloc1(numLeaves, &ranksNew));
46 PetscCall(PetscMalloc1(numLeaves, &localPointsNew));
47 PetscCall(PetscMalloc1(numLeaves, &remotePointsNew));
48 for (l = 0; l < numLeaves; ++l) {
49 ranksNew[l] = ranks[l];
50 localPointsNew[l] = l;
51 remotePointsNew[l].index = 0;
52 remotePointsNew[l].rank = ranksNew[l];
53 }
54 PetscCall(PetscFree(ranks));
55 if (processRanks) PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)dm), numLeaves, ranksNew, PETSC_OWN_POINTER, processRanks));
56 else PetscCall(PetscFree(ranksNew));
57 if (sfProcess) {
58 PetscCall(PetscSFCreate(PetscObjectComm((PetscObject)dm), sfProcess));
59 PetscCall(PetscObjectSetName((PetscObject)*sfProcess, "Process SF"));
60 PetscCall(PetscSFSetFromOptions(*sfProcess));
61 PetscCall(PetscSFSetGraph(*sfProcess, size, numLeaves, localPointsNew, PETSC_OWN_POINTER, remotePointsNew, PETSC_OWN_POINTER));
62 }
63 PetscFunctionReturn(PETSC_SUCCESS);
64 }
65
66 /*@
67 DMPlexCreateCoarsePointIS - Creates an `IS` covering the coarse `DM` chart with the fine points as data
68
69 Collective
70
71 Input Parameter:
72 . dm - The coarse `DM`
73
74 Output Parameter:
75 . fpointIS - The `IS` of all the fine points which exist in the original coarse mesh
76
77 Level: developer
78
79 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `IS`, `DMRefine()`, `DMPlexSetRefinementUniform()`, `DMPlexGetSubpointIS()`
80 @*/
DMPlexCreateCoarsePointIS(DM dm,IS * fpointIS)81 PetscErrorCode DMPlexCreateCoarsePointIS(DM dm, IS *fpointIS)
82 {
83 DMPlexTransform tr;
84 PetscInt *fpoints;
85 PetscInt pStart, pEnd, p, vStart, vEnd, v;
86
87 PetscFunctionBegin;
88 PetscCall(DMPlexGetChart(dm, &pStart, &pEnd));
89 PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd));
90 PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), &tr));
91 PetscCall(DMPlexTransformSetUp(tr));
92 PetscCall(PetscMalloc1(pEnd - pStart, &fpoints));
93 for (p = 0; p < pEnd - pStart; ++p) fpoints[p] = -1;
94 for (v = vStart; v < vEnd; ++v) {
95 PetscInt vNew = -1; /* quiet overzealous may be used uninitialized check */
96
97 PetscCall(DMPlexTransformGetTargetPoint(tr, DM_POLYTOPE_POINT, DM_POLYTOPE_POINT, p, 0, &vNew));
98 fpoints[v - pStart] = vNew;
99 }
100 PetscCall(DMPlexTransformDestroy(&tr));
101 PetscCall(ISCreateGeneral(PETSC_COMM_SELF, pEnd - pStart, fpoints, PETSC_OWN_POINTER, fpointIS));
102 PetscFunctionReturn(PETSC_SUCCESS);
103 }
104
105 /*@
106 DMPlexSetTransformType - Set the transform type for uniform refinement
107
108 Input Parameters:
109 + dm - The `DM`
110 - type - The transform type for uniform refinement
111
112 Level: developer
113
114 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexTransformType`, `DMRefine()`, `DMPlexGetTransformType()`, `DMPlexSetRefinementUniform()`
115 @*/
DMPlexSetTransformType(DM dm,DMPlexTransformType type)116 PetscErrorCode DMPlexSetTransformType(DM dm, DMPlexTransformType type)
117 {
118 DM_Plex *mesh = (DM_Plex *)dm->data;
119
120 PetscFunctionBegin;
121 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMPLEX);
122 if (type) PetscAssertPointer(type, 2);
123 PetscCall(PetscFree(mesh->transformType));
124 PetscCall(PetscStrallocpy(type, &mesh->transformType));
125 PetscFunctionReturn(PETSC_SUCCESS);
126 }
127
128 /*@C
129 DMPlexGetTransformType - Retrieve the transform type for uniform refinement
130
131 Input Parameter:
132 . dm - The `DM`
133
134 Output Parameter:
135 . type - The transform type for uniform refinement
136
137 Level: developer
138
139 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexTransformType`, `DMRefine()`, `DMPlexSetTransformType()`, `DMPlexGetRefinementUniform()`
140 @*/
DMPlexGetTransformType(DM dm,DMPlexTransformType * type)141 PetscErrorCode DMPlexGetTransformType(DM dm, DMPlexTransformType *type)
142 {
143 DM_Plex *mesh = (DM_Plex *)dm->data;
144
145 PetscFunctionBegin;
146 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMPLEX);
147 PetscAssertPointer(type, 2);
148 *type = mesh->transformType;
149 PetscFunctionReturn(PETSC_SUCCESS);
150 }
151
DMPlexSetTransform(DM dm,DMPlexTransform tr)152 PetscErrorCode DMPlexSetTransform(DM dm, DMPlexTransform tr)
153 {
154 DM_Plex *mesh = (DM_Plex *)dm->data;
155
156 PetscFunctionBegin;
157 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMPLEX);
158 if (tr) PetscValidHeaderSpecific(tr, DMPLEXTRANSFORM_CLASSID, 2);
159 PetscCall(PetscObjectReference((PetscObject)tr));
160 PetscCall(DMPlexTransformDestroy(&mesh->transform));
161 mesh->transform = tr;
162 PetscFunctionReturn(PETSC_SUCCESS);
163 }
164
DMPlexGetTransform(DM dm,DMPlexTransform * tr)165 PetscErrorCode DMPlexGetTransform(DM dm, DMPlexTransform *tr)
166 {
167 DM_Plex *mesh = (DM_Plex *)dm->data;
168
169 PetscFunctionBegin;
170 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMPLEX);
171 PetscAssertPointer(tr, 2);
172 *tr = mesh->transform;
173 PetscFunctionReturn(PETSC_SUCCESS);
174 }
175
DMPlexSetSaveTransform(DM dm,PetscBool save)176 PetscErrorCode DMPlexSetSaveTransform(DM dm, PetscBool save)
177 {
178 DM_Plex *mesh = (DM_Plex *)dm->data;
179
180 PetscFunctionBegin;
181 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMPLEX);
182 mesh->saveTransform = save;
183 PetscFunctionReturn(PETSC_SUCCESS);
184 }
185
DMPlexGetSaveTransform(DM dm,PetscBool * save)186 PetscErrorCode DMPlexGetSaveTransform(DM dm, PetscBool *save)
187 {
188 DM_Plex *mesh = (DM_Plex *)dm->data;
189
190 PetscFunctionBegin;
191 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMPLEX);
192 PetscAssertPointer(save, 2);
193 *save = mesh->saveTransform;
194 PetscFunctionReturn(PETSC_SUCCESS);
195 }
196
197 /*@
198 DMPlexSetRefinementUniform - Set the flag for uniform refinement
199
200 Input Parameters:
201 + dm - The `DM`
202 - refinementUniform - The flag for uniform refinement
203
204 Level: developer
205
206 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMRefine()`, `DMPlexGetRefinementUniform()`, `DMPlexGetRefinementLimit()`, `DMPlexSetRefinementLimit()`
207 @*/
DMPlexSetRefinementUniform(DM dm,PetscBool refinementUniform)208 PetscErrorCode DMPlexSetRefinementUniform(DM dm, PetscBool refinementUniform)
209 {
210 DM_Plex *mesh = (DM_Plex *)dm->data;
211
212 PetscFunctionBegin;
213 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMPLEX);
214 mesh->refinementUniform = refinementUniform;
215 PetscFunctionReturn(PETSC_SUCCESS);
216 }
217
218 /*@
219 DMPlexGetRefinementUniform - Retrieve the flag for uniform refinement
220
221 Input Parameter:
222 . dm - The `DM`
223
224 Output Parameter:
225 . refinementUniform - The flag for uniform refinement
226
227 Level: developer
228
229 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMRefine()`, `DMPlexSetRefinementUniform()`, `DMPlexGetRefinementLimit()`, `DMPlexSetRefinementLimit()`
230 @*/
DMPlexGetRefinementUniform(DM dm,PetscBool * refinementUniform)231 PetscErrorCode DMPlexGetRefinementUniform(DM dm, PetscBool *refinementUniform)
232 {
233 DM_Plex *mesh = (DM_Plex *)dm->data;
234
235 PetscFunctionBegin;
236 PetscValidHeaderSpecificType(dm, DM_CLASSID, 1, DMPLEX);
237 PetscAssertPointer(refinementUniform, 2);
238 *refinementUniform = mesh->refinementUniform;
239 PetscFunctionReturn(PETSC_SUCCESS);
240 }
241
242 /*@
243 DMPlexSetRefinementLimit - Set the maximum cell volume for refinement
244
245 Input Parameters:
246 + dm - The `DM`
247 - refinementLimit - The maximum cell volume in the refined mesh
248
249 Level: developer
250
251 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMRefine()`, `DMPlexGetRefinementLimit()`, `DMPlexGetRefinementUniform()`, `DMPlexSetRefinementUniform()`
252 @*/
DMPlexSetRefinementLimit(DM dm,PetscReal refinementLimit)253 PetscErrorCode DMPlexSetRefinementLimit(DM dm, PetscReal refinementLimit)
254 {
255 DM_Plex *mesh = (DM_Plex *)dm->data;
256
257 PetscFunctionBegin;
258 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
259 mesh->refinementLimit = refinementLimit;
260 PetscFunctionReturn(PETSC_SUCCESS);
261 }
262
263 /*@
264 DMPlexGetRefinementLimit - Retrieve the maximum cell volume for refinement
265
266 Input Parameter:
267 . dm - The `DM`
268
269 Output Parameter:
270 . refinementLimit - The maximum cell volume in the refined mesh
271
272 Level: developer
273
274 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMRefine()`, `DMPlexSetRefinementLimit()`, `DMPlexGetRefinementUniform()`, `DMPlexSetRefinementUniform()`
275 @*/
DMPlexGetRefinementLimit(DM dm,PetscReal * refinementLimit)276 PetscErrorCode DMPlexGetRefinementLimit(DM dm, PetscReal *refinementLimit)
277 {
278 DM_Plex *mesh = (DM_Plex *)dm->data;
279
280 PetscFunctionBegin;
281 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
282 PetscAssertPointer(refinementLimit, 2);
283 /* if (mesh->refinementLimit < 0) = getMaxVolume()/2.0; */
284 *refinementLimit = mesh->refinementLimit;
285 PetscFunctionReturn(PETSC_SUCCESS);
286 }
287
288 /*@C
289 DMPlexSetRefinementFunction - Set the function giving the maximum cell volume for refinement
290
291 Input Parameters:
292 + dm - The `DM`
293 - refinementFunc - Function giving the maximum cell volume in the refined mesh
294
295 Calling Sequence of `refinementFunc`:
296 + coords - Coordinates of the current point, usually a cell centroid
297 - limit - The maximum cell volume for a cell containing this point
298
299 Level: developer
300
301 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMRefine()`, `DMPlexGetRefinementFunction()`, `DMPlexGetRefinementUniform()`, `DMPlexSetRefinementUniform()`, `DMPlexGetRefinementLimit()`, `DMPlexSetRefinementLimit()`
302 @*/
DMPlexSetRefinementFunction(DM dm,PetscErrorCode (* refinementFunc)(const PetscReal coords[],PetscReal * limit))303 PetscErrorCode DMPlexSetRefinementFunction(DM dm, PetscErrorCode (*refinementFunc)(const PetscReal coords[], PetscReal *limit))
304 {
305 DM_Plex *mesh = (DM_Plex *)dm->data;
306
307 PetscFunctionBegin;
308 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
309 mesh->refinementFunc = refinementFunc;
310 PetscFunctionReturn(PETSC_SUCCESS);
311 }
312
313 /*@C
314 DMPlexGetRefinementFunction - Get the function giving the maximum cell volume for refinement
315
316 Input Parameter:
317 . dm - The `DM`
318
319 Output Parameter:
320 . refinementFunc - Function giving the maximum cell volume in the refined mesh
321
322 Calling Sequence of `refinementFunc`:
323 + coords - Coordinates of the current point, usually a cell centroid
324 - limit - The maximum cell volume for a cell containing this point
325
326 Level: developer
327
328 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMRefine()`, `DMPlexSetRefinementFunction()`, `DMPlexGetRefinementUniform()`, `DMPlexSetRefinementUniform()`, `DMPlexGetRefinementLimit()`, `DMPlexSetRefinementLimit()`
329 @*/
DMPlexGetRefinementFunction(DM dm,PetscErrorCode (** refinementFunc)(const PetscReal coords[],PetscReal * limit))330 PetscErrorCode DMPlexGetRefinementFunction(DM dm, PetscErrorCode (**refinementFunc)(const PetscReal coords[], PetscReal *limit))
331 {
332 DM_Plex *mesh = (DM_Plex *)dm->data;
333
334 PetscFunctionBegin;
335 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
336 PetscAssertPointer(refinementFunc, 2);
337 *refinementFunc = mesh->refinementFunc;
338 PetscFunctionReturn(PETSC_SUCCESS);
339 }
340
DMRefine_Plex(DM dm,MPI_Comm comm,DM * rdm)341 PetscErrorCode DMRefine_Plex(DM dm, MPI_Comm comm, DM *rdm)
342 {
343 PetscBool isUniform;
344
345 PetscFunctionBegin;
346 PetscCall(DMPlexGetRefinementUniform(dm, &isUniform));
347 PetscCall(DMViewFromOptions(dm, NULL, "-initref_dm_view"));
348 if (isUniform) {
349 DMPlexTransform tr;
350 DM cdm, rcdm;
351 DMPlexTransformType trType;
352 const char *prefix;
353 PetscOptions options;
354 PetscInt cDegree;
355 PetscBool useCeed, save;
356
357 PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), &tr));
358 PetscCall(DMPlexTransformSetDM(tr, dm));
359 PetscCall(DMPlexGetTransformType(dm, &trType));
360 if (trType) PetscCall(DMPlexTransformSetType(tr, trType));
361 PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix));
362 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)tr, prefix));
363 PetscCall(PetscObjectGetOptions((PetscObject)dm, &options));
364 PetscCall(PetscObjectSetOptions((PetscObject)tr, options));
365 PetscCall(DMPlexTransformSetFromOptions(tr));
366 PetscCall(PetscObjectSetOptions((PetscObject)tr, NULL));
367 PetscCall(DMPlexTransformSetUp(tr));
368 PetscCall(PetscObjectViewFromOptions((PetscObject)tr, NULL, "-dm_plex_transform_view"));
369 PetscCall(DMPlexTransformApply(tr, dm, rdm));
370 PetscCall(DMPlexSetRegularRefinement(*rdm, PETSC_TRUE));
371 PetscCall(DMPlexGetUseCeed(dm, &useCeed));
372 PetscCall(DMPlexSetUseCeed(*rdm, useCeed));
373 PetscCall(DMSetMatType(*rdm, dm->mattype));
374 PetscCall(DMCopyDisc(dm, *rdm));
375 PetscCall(DMGetCoordinateDM(dm, &cdm));
376 PetscCall(DMGetCoordinateDM(*rdm, &rcdm));
377 PetscCall(DMGetCoordinateDegree_Internal(dm, &cDegree));
378 {
379 PetscDS cds, rcds;
380
381 PetscCall(DMPlexCreateCoordinateSpace(*rdm, cDegree, PETSC_FALSE, PETSC_TRUE));
382 PetscCall(DMGetCoordinateDM(*rdm, &rcdm));
383 PetscCall(DMGetDS(cdm, &cds));
384 PetscCall(DMGetDS(rcdm, &rcds));
385 PetscCall(PetscDSCopyConstants(cds, rcds));
386 }
387 PetscCall(DMPlexGetUseCeed(cdm, &useCeed));
388 PetscCall(DMPlexSetUseCeed(rcdm, useCeed));
389 if (useCeed) {
390 PetscCall(DMPlexSetUseMatClosurePermutation(rcdm, PETSC_FALSE));
391 PetscCall(DMUseTensorOrder(rcdm, PETSC_TRUE));
392 }
393 PetscCall(DMPlexTransformCreateDiscLabels(tr, *rdm));
394 PetscCall(DMPlexGetSaveTransform(dm, &save));
395 if (save) PetscCall(DMPlexSetTransform(*rdm, tr));
396 PetscCall(DMPlexTransformDestroy(&tr));
397 } else {
398 PetscCall(DMPlexRefine_Internal(dm, NULL, NULL, NULL, rdm));
399 }
400 if (*rdm) {
401 ((DM_Plex *)(*rdm)->data)->printFEM = ((DM_Plex *)dm->data)->printFEM;
402 ((DM_Plex *)(*rdm)->data)->printL2 = ((DM_Plex *)dm->data)->printL2;
403 }
404 PetscCall(DMViewFromOptions(*rdm, NULL, "-postref_dm_view"));
405 PetscFunctionReturn(PETSC_SUCCESS);
406 }
407
DMRefineHierarchy_Plex(DM dm,PetscInt nlevels,DM rdm[])408 PetscErrorCode DMRefineHierarchy_Plex(DM dm, PetscInt nlevels, DM rdm[])
409 {
410 DM cdm = dm;
411 PetscInt r;
412 PetscBool isUniform, localized, useCeed;
413
414 PetscFunctionBegin;
415 PetscCall(DMPlexGetRefinementUniform(dm, &isUniform));
416 PetscCall(DMGetCoordinatesLocalized(dm, &localized));
417 if (isUniform) {
418 for (r = 0; r < nlevels; ++r) {
419 DMPlexTransform tr;
420 DM codm, rcodm;
421 const char *prefix;
422
423 PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)cdm), &tr));
424 PetscCall(PetscObjectGetOptionsPrefix((PetscObject)cdm, &prefix));
425 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)tr, prefix));
426 PetscCall(DMPlexTransformSetDM(tr, cdm));
427 PetscCall(DMPlexTransformSetFromOptions(tr));
428 PetscCall(DMPlexTransformSetUp(tr));
429 PetscCall(DMPlexTransformApply(tr, cdm, &rdm[r]));
430 PetscCall(DMSetCoarsenLevel(rdm[r], cdm->leveldown));
431 PetscCall(DMSetRefineLevel(rdm[r], cdm->levelup + 1));
432 PetscCall(DMSetMatType(rdm[r], dm->mattype));
433 PetscCall(DMPlexGetUseCeed(dm, &useCeed));
434 PetscCall(DMPlexSetUseCeed(rdm[r], useCeed));
435 PetscCall(DMCopyDisc(cdm, rdm[r]));
436 PetscCall(DMGetCoordinateDM(dm, &codm));
437 PetscCall(DMGetCoordinateDM(rdm[r], &rcodm));
438 PetscCall(DMCopyDisc(codm, rcodm));
439 PetscCall(DMPlexGetUseCeed(codm, &useCeed));
440 PetscCall(DMPlexSetUseCeed(rcodm, useCeed));
441 if (useCeed) {
442 PetscCall(DMPlexSetUseMatClosurePermutation(rcodm, PETSC_FALSE));
443 PetscCall(DMUseTensorOrder(rcodm, PETSC_TRUE));
444 }
445 PetscCall(DMPlexTransformCreateDiscLabels(tr, rdm[r]));
446 PetscCall(DMSetCoarseDM(rdm[r], cdm));
447 PetscCall(DMPlexSetRegularRefinement(rdm[r], PETSC_TRUE));
448 if (rdm[r]) {
449 ((DM_Plex *)rdm[r]->data)->printFEM = ((DM_Plex *)dm->data)->printFEM;
450 ((DM_Plex *)rdm[r]->data)->printL2 = ((DM_Plex *)dm->data)->printL2;
451 }
452 cdm = rdm[r];
453 PetscCall(DMPlexTransformDestroy(&tr));
454 }
455 } else {
456 for (r = 0; r < nlevels; ++r) {
457 PetscCall(DMRefine(cdm, PetscObjectComm((PetscObject)dm), &rdm[r]));
458 PetscCall(DMPlexGetUseCeed(dm, &useCeed));
459 PetscCall(DMPlexSetUseCeed(rdm[r], useCeed));
460 PetscCall(DMCopyDisc(cdm, rdm[r]));
461 if (localized) PetscCall(DMLocalizeCoordinates(rdm[r]));
462 PetscCall(DMSetCoarseDM(rdm[r], cdm));
463 if (rdm[r]) {
464 ((DM_Plex *)rdm[r]->data)->printFEM = ((DM_Plex *)dm->data)->printFEM;
465 ((DM_Plex *)rdm[r]->data)->printL2 = ((DM_Plex *)dm->data)->printL2;
466 }
467 cdm = rdm[r];
468 }
469 }
470 PetscFunctionReturn(PETSC_SUCCESS);
471 }
472