Actual source code: forest.c

  1: #include <petsc/private/dmforestimpl.h>
  2: #include <petsc/private/dmimpl.h>
  3: #include <petsc/private/dmlabelimpl.h>
  4: #include <petscsf.h>

  6: PetscBool DMForestPackageInitialized = PETSC_FALSE;

  8: typedef struct _DMForestTypeLink *DMForestTypeLink;

 10: struct _DMForestTypeLink {
 11:   char            *name;
 12:   DMForestTypeLink next;
 13: };

 15: DMForestTypeLink DMForestTypeList;

 17: static PetscErrorCode DMForestPackageFinalize(void)
 18: {
 19:   DMForestTypeLink oldLink, link = DMForestTypeList;

 21:   PetscFunctionBegin;
 22:   while (link) {
 23:     oldLink = link;
 24:     PetscCall(PetscFree(oldLink->name));
 25:     link = oldLink->next;
 26:     PetscCall(PetscFree(oldLink));
 27:   }
 28:   PetscFunctionReturn(PETSC_SUCCESS);
 29: }

 31: static PetscErrorCode DMForestPackageInitialize(void)
 32: {
 33:   PetscFunctionBegin;
 34:   if (DMForestPackageInitialized) PetscFunctionReturn(PETSC_SUCCESS);
 35:   DMForestPackageInitialized = PETSC_TRUE;

 37:   PetscCall(DMForestRegisterType(DMFOREST));
 38:   PetscCall(PetscRegisterFinalize(DMForestPackageFinalize));
 39:   PetscFunctionReturn(PETSC_SUCCESS);
 40: }

 42: /*@C
 43:   DMForestRegisterType - Registers a `DMType` as a subtype of `DMFOREST` (so that `DMIsForest()` will be correct)

 45:   Not Collective

 47:   Input Parameter:
 48: . name - the name of the type

 50:   Level: advanced

 52: .seealso: `DMFOREST`, `DMIsForest()`
 53: @*/
 54: PetscErrorCode DMForestRegisterType(DMType name)
 55: {
 56:   DMForestTypeLink link;

 58:   PetscFunctionBegin;
 59:   PetscCall(DMForestPackageInitialize());
 60:   PetscCall(PetscNew(&link));
 61:   PetscCall(PetscStrallocpy(name, &link->name));
 62:   link->next       = DMForestTypeList;
 63:   DMForestTypeList = link;
 64:   PetscFunctionReturn(PETSC_SUCCESS);
 65: }

 67: /*@
 68:   DMIsForest - Check whether a DM uses the DMFOREST interface for hierarchically-refined meshes

 70:   Not Collective

 72:   Input Parameter:
 73: . dm - the DM object

 75:   Output Parameter:
 76: . isForest - whether dm is a subtype of DMFOREST

 78:   Level: intermediate

 80: .seealso: `DMFOREST`, `DMForestRegisterType()`
 81: @*/
 82: PetscErrorCode DMIsForest(DM dm, PetscBool *isForest)
 83: {
 84:   DMForestTypeLink link = DMForestTypeList;

 86:   PetscFunctionBegin;
 87:   while (link) {
 88:     PetscBool sameType;
 89:     PetscCall(PetscObjectTypeCompare((PetscObject)dm, link->name, &sameType));
 90:     if (sameType) {
 91:       *isForest = PETSC_TRUE;
 92:       PetscFunctionReturn(PETSC_SUCCESS);
 93:     }
 94:     link = link->next;
 95:   }
 96:   *isForest = PETSC_FALSE;
 97:   PetscFunctionReturn(PETSC_SUCCESS);
 98: }

100: /*@
101:   DMForestTemplate - Create a new `DM` that will be adapted from a source `DM`.

103:   Collective

105:   Input Parameters:
106: + dm   - the source `DM` object
107: - comm - the communicator for the new `DM` (this communicator is currently ignored, but is present so that `DMForestTemplate()` can be used within `DMCoarsen()`)

109:   Output Parameter:
110: . tdm - the new `DM` object

112:   Level: intermediate

114:   Notes:
115:   The new `DM` reproduces the configuration of the source, but is not yet setup, so that the
116:   user can then define only the ways that the new `DM` should differ (by, e.g., refinement or
117:   repartitioning).  The source `DM` is also set as the adaptivity source `DM` of the new `DM`
118:   (see `DMForestSetAdaptivityForest()`).

120: .seealso: `DM`, `DMFOREST`, `DMForestSetAdaptivityForest()`
121: @*/
122: PetscErrorCode DMForestTemplate(DM dm, MPI_Comm comm, DM *tdm)
123: {
124:   DM_Forest                 *forest = (DM_Forest *)dm->data;
125:   DMType                     type;
126:   DM                         base;
127:   DMForestTopology           topology;
128:   MatType                    mtype;
129:   PetscInt                   dim, overlap, ref, factor;
130:   DMForestAdaptivityStrategy strat;
131:   void                      *ctx;
132:   PetscErrorCode (*map)(DM, PetscInt, PetscInt, const PetscReal[], PetscReal[], void *);
133:   void *mapCtx;

135:   PetscFunctionBegin;
137:   PetscCall(DMCreate(PetscObjectComm((PetscObject)dm), tdm));
138:   PetscCall(DMGetType(dm, &type));
139:   PetscCall(DMSetType(*tdm, type));
140:   PetscCall(DMForestGetBaseDM(dm, &base));
141:   PetscCall(DMForestSetBaseDM(*tdm, base));
142:   PetscCall(DMForestGetTopology(dm, &topology));
143:   PetscCall(DMForestSetTopology(*tdm, topology));
144:   PetscCall(DMForestGetAdjacencyDimension(dm, &dim));
145:   PetscCall(DMForestSetAdjacencyDimension(*tdm, dim));
146:   PetscCall(DMForestGetPartitionOverlap(dm, &overlap));
147:   PetscCall(DMForestSetPartitionOverlap(*tdm, overlap));
148:   PetscCall(DMForestGetMinimumRefinement(dm, &ref));
149:   PetscCall(DMForestSetMinimumRefinement(*tdm, ref));
150:   PetscCall(DMForestGetMaximumRefinement(dm, &ref));
151:   PetscCall(DMForestSetMaximumRefinement(*tdm, ref));
152:   PetscCall(DMForestGetAdaptivityStrategy(dm, &strat));
153:   PetscCall(DMForestSetAdaptivityStrategy(*tdm, strat));
154:   PetscCall(DMForestGetGradeFactor(dm, &factor));
155:   PetscCall(DMForestSetGradeFactor(*tdm, factor));
156:   PetscCall(DMForestGetBaseCoordinateMapping(dm, &map, &mapCtx));
157:   PetscCall(DMForestSetBaseCoordinateMapping(*tdm, map, mapCtx));
158:   if (forest->ftemplate) PetscCall((*forest->ftemplate)(dm, *tdm));
159:   PetscCall(DMForestSetAdaptivityForest(*tdm, dm));
160:   PetscCall(DMCopyDisc(dm, *tdm));
161:   PetscCall(DMGetApplicationContext(dm, &ctx));
162:   PetscCall(DMSetApplicationContext(*tdm, &ctx));
163:   {
164:     const PetscReal *maxCell, *L, *Lstart;

166:     PetscCall(DMGetPeriodicity(dm, &maxCell, &Lstart, &L));
167:     PetscCall(DMSetPeriodicity(*tdm, maxCell, Lstart, L));
168:   }
169:   PetscCall(DMGetMatType(dm, &mtype));
170:   PetscCall(DMSetMatType(*tdm, mtype));
171:   PetscFunctionReturn(PETSC_SUCCESS);
172: }

174: static PetscErrorCode DMInitialize_Forest(DM dm);

176: PETSC_EXTERN PetscErrorCode DMClone_Forest(DM dm, DM *newdm)
177: {
178:   DM_Forest  *forest = (DM_Forest *)dm->data;
179:   const char *type;

181:   PetscFunctionBegin;
182:   forest->refct++;
183:   (*newdm)->data = forest;
184:   PetscCall(PetscObjectGetType((PetscObject)dm, &type));
185:   PetscCall(PetscObjectChangeTypeName((PetscObject)*newdm, type));
186:   PetscCall(DMInitialize_Forest(*newdm));
187:   PetscFunctionReturn(PETSC_SUCCESS);
188: }

190: static PetscErrorCode DMDestroy_Forest(DM dm)
191: {
192:   DM_Forest *forest = (DM_Forest *)dm->data;

194:   PetscFunctionBegin;
195:   PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMConvert_plex_p4est_C", NULL));
196:   PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMConvert_p4est_plex_C", NULL));
197:   PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMConvert_plex_p8est_C", NULL));
198:   PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMConvert_p8est_plex_C", NULL));
199:   PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMSetUpGLVisViewer_C", NULL));
200:   PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMCreateNeumannOverlap_C", NULL));
201:   PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMPlexGetOverlap_C", NULL));
202:   PetscCall(PetscObjectComposeFunction((PetscObject)dm, "MatComputeNeumannOverlap_C", NULL));
203:   if (--forest->refct > 0) PetscFunctionReturn(PETSC_SUCCESS);
204:   if (forest->destroy) PetscCall((*forest->destroy)(dm));
205:   PetscCall(PetscSFDestroy(&forest->cellSF));
206:   PetscCall(PetscSFDestroy(&forest->preCoarseToFine));
207:   PetscCall(PetscSFDestroy(&forest->coarseToPreFine));
208:   PetscCall(DMLabelDestroy(&forest->adaptLabel));
209:   PetscCall(PetscFree(forest->adaptStrategy));
210:   PetscCall(DMDestroy(&forest->base));
211:   PetscCall(DMDestroy(&forest->adapt));
212:   PetscCall(PetscFree(forest->topology));
213:   PetscCall(PetscFree(forest));
214:   PetscFunctionReturn(PETSC_SUCCESS);
215: }

217: /*@
218:   DMForestSetTopology - Set the topology of a `DMFOREST` during the pre-setup phase.  The topology is a string (e.g.
219:   "cube", "shell") and can be interpreted by subtypes of `DMFOREST`) to construct the base DM of a forest during
220:   `DMSetUp()`.

222:   Logically collectiv

224:   Input Parameters:
225: + dm       - the forest
226: - topology - the topology of the forest

228:   Level: intermediate

230: .seealso: `DM`, `DMFOREST`, `DMForestGetTopology()`, `DMForestSetBaseDM()`
231: @*/
232: PetscErrorCode DMForestSetTopology(DM dm, DMForestTopology topology)
233: {
234:   DM_Forest *forest = (DM_Forest *)dm->data;

236:   PetscFunctionBegin;
238:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the topology after setup");
239:   PetscCall(PetscFree(forest->topology));
240:   PetscCall(PetscStrallocpy((const char *)topology, (char **)&forest->topology));
241:   PetscFunctionReturn(PETSC_SUCCESS);
242: }

244: /*@
245:   DMForestGetTopology - Get a string describing the topology of a `DMFOREST`.

247:   Not Collective

249:   Input Parameter:
250: . dm - the forest

252:   Output Parameter:
253: . topology - the topology of the forest (e.g., 'cube', 'shell')

255:   Level: intermediate

257: .seealso: `DM`, `DMFOREST`, `DMForestSetTopology()`
258: @*/
259: PetscErrorCode DMForestGetTopology(DM dm, DMForestTopology *topology)
260: {
261:   DM_Forest *forest = (DM_Forest *)dm->data;

263:   PetscFunctionBegin;
265:   PetscAssertPointer(topology, 2);
266:   *topology = forest->topology;
267:   PetscFunctionReturn(PETSC_SUCCESS);
268: }

270: /*@
271:   DMForestSetBaseDM - During the pre-setup phase, set the `DM` that defines the base mesh of a
272:   `DMFOREST` forest.

274:   Logically Collective

276:   Input Parameters:
277: + dm   - the forest
278: - base - the base `DM` of the forest

280:   Level: intermediate

282:   Notes:
283:   The forest will be hierarchically refined from the base, and all refinements/coarsenings of
284:   the forest will share its base.  In general, two forest must share a base to be comparable,
285:   to do things like construct interpolators.

287:   Currently the base `DM` must be a `DMPLEX`

289: .seealso: `DM`, `DMFOREST`, `DMForestGetBaseDM()`
290: @*/
291: PetscErrorCode DMForestSetBaseDM(DM dm, DM base)
292: {
293:   DM_Forest *forest = (DM_Forest *)dm->data;
294:   PetscInt   dim, dimEmbed;

296:   PetscFunctionBegin;
298:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the base after setup");
299:   PetscCall(PetscObjectReference((PetscObject)base));
300:   PetscCall(DMDestroy(&forest->base));
301:   forest->base = base;
302:   if (base) {
303:     const PetscReal *maxCell, *Lstart, *L;

306:     PetscCall(DMGetDimension(base, &dim));
307:     PetscCall(DMSetDimension(dm, dim));
308:     PetscCall(DMGetCoordinateDim(base, &dimEmbed));
309:     PetscCall(DMSetCoordinateDim(dm, dimEmbed));
310:     PetscCall(DMGetPeriodicity(base, &maxCell, &Lstart, &L));
311:     PetscCall(DMSetPeriodicity(dm, maxCell, Lstart, L));
312:   } else PetscCall(DMSetPeriodicity(dm, NULL, NULL, NULL));
313:   PetscFunctionReturn(PETSC_SUCCESS);
314: }

316: /*@
317:   DMForestGetBaseDM - Get the base `DM` of a `DMFOREST`

319:   Not Collective

321:   Input Parameter:
322: . dm - the forest

324:   Output Parameter:
325: . base - the base `DM` of the forest

327:   Level: intermediate

329:   Notes:
330:   After DMSetUp(), the base DM will be redundantly distributed across MPI processes

332:   The forest will be hierarchically refined from the base, and all refinements/coarsenings of
333:   the forest will share its base.  In general, two forest must share a base to be comparable,
334:   to do things like construct interpolators.

336: .seealso: `DM`, `DMFOREST`, `DMForestSetBaseDM()`
337: @*/
338: PetscErrorCode DMForestGetBaseDM(DM dm, DM *base)
339: {
340:   DM_Forest *forest = (DM_Forest *)dm->data;

342:   PetscFunctionBegin;
344:   PetscAssertPointer(base, 2);
345:   *base = forest->base;
346:   PetscFunctionReturn(PETSC_SUCCESS);
347: }

349: PetscErrorCode DMForestSetBaseCoordinateMapping(DM dm, PetscErrorCode (*func)(DM, PetscInt, PetscInt, const PetscReal[], PetscReal[], void *), void *ctx)
350: {
351:   DM_Forest *forest = (DM_Forest *)dm->data;

353:   PetscFunctionBegin;
355:   forest->mapcoordinates    = func;
356:   forest->mapcoordinatesctx = ctx;
357:   PetscFunctionReturn(PETSC_SUCCESS);
358: }

360: PetscErrorCode DMForestGetBaseCoordinateMapping(DM dm, PetscErrorCode (**func)(DM, PetscInt, PetscInt, const PetscReal[], PetscReal[], void *), void *ctx)
361: {
362:   DM_Forest *forest = (DM_Forest *)dm->data;

364:   PetscFunctionBegin;
366:   if (func) *func = forest->mapcoordinates;
367:   if (ctx) *((void **)ctx) = forest->mapcoordinatesctx;
368:   PetscFunctionReturn(PETSC_SUCCESS);
369: }

371: /*@
372:   DMForestSetAdaptivityForest - During the pre-setup phase, set the forest from which the
373:   current forest will be adapted (e.g., the current forest will be
374:   refined/coarsened/repartitioned from it) in `DMSetUp()`.

376:   Logically Collective

378:   Input Parameters:
379: + dm    - the new forest, which will be constructed from adapt
380: - adapt - the old forest

382:   Level: intermediate

384:   Note:
385:   Usually not needed by users directly, `DMForestTemplate()` constructs a new forest to be
386:   adapted from an old forest and calls this routine.

388:   This can be called after setup with `adapt` = `NULL`, which will clear all internal data
389:   related to the adaptivity forest from `dm`. This way, repeatedly adapting does not leave
390:   stale `DM` objects in memory.

392: .seealso: `DM`, `DMFOREST`, `DMForestGetAdaptivityForest()`, `DMForestSetAdaptivityPurpose()`
393: @*/
394: PetscErrorCode DMForestSetAdaptivityForest(DM dm, DM adapt)
395: {
396:   DM_Forest *forest, *adaptForest, *oldAdaptForest;
397:   DM         oldAdapt;
398:   PetscBool  isForest;

400:   PetscFunctionBegin;
403:   PetscCall(DMIsForest(dm, &isForest));
404:   if (!isForest) PetscFunctionReturn(PETSC_SUCCESS);
405:   PetscCheck(adapt == NULL || !dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the adaptation forest after setup");
406:   forest = (DM_Forest *)dm->data;
407:   PetscCall(DMForestGetAdaptivityForest(dm, &oldAdapt));
408:   adaptForest    = (DM_Forest *)(adapt ? adapt->data : NULL);
409:   oldAdaptForest = (DM_Forest *)(oldAdapt ? oldAdapt->data : NULL);
410:   if (adaptForest != oldAdaptForest) {
411:     PetscCall(PetscSFDestroy(&forest->preCoarseToFine));
412:     PetscCall(PetscSFDestroy(&forest->coarseToPreFine));
413:     if (forest->clearadaptivityforest) PetscCall((*forest->clearadaptivityforest)(dm));
414:   }
415:   switch (forest->adaptPurpose) {
416:   case DM_ADAPT_DETERMINE:
417:     PetscCall(PetscObjectReference((PetscObject)adapt));
418:     PetscCall(DMDestroy(&forest->adapt));
419:     forest->adapt = adapt;
420:     break;
421:   case DM_ADAPT_REFINE:
422:     PetscCall(DMSetCoarseDM(dm, adapt));
423:     break;
424:   case DM_ADAPT_COARSEN:
425:   case DM_ADAPT_COARSEN_LAST:
426:     PetscCall(DMSetFineDM(dm, adapt));
427:     break;
428:   default:
429:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "invalid adaptivity purpose");
430:   }
431:   PetscFunctionReturn(PETSC_SUCCESS);
432: }

434: /*@
435:   DMForestGetAdaptivityForest - Get the forest from which the current forest is adapted.

437:   Not Collective

439:   Input Parameter:
440: . dm - the forest

442:   Output Parameter:
443: . adapt - the forest from which `dm` is/was adapted

445:   Level: intermediate

447: .seealso: `DM`, `DMFOREST`, `DMForestSetAdaptivityForest()`, `DMForestSetAdaptivityPurpose()`
448: @*/
449: PetscErrorCode DMForestGetAdaptivityForest(DM dm, DM *adapt)
450: {
451:   DM_Forest *forest;

453:   PetscFunctionBegin;
455:   forest = (DM_Forest *)dm->data;
456:   switch (forest->adaptPurpose) {
457:   case DM_ADAPT_DETERMINE:
458:     *adapt = forest->adapt;
459:     break;
460:   case DM_ADAPT_REFINE:
461:     PetscCall(DMGetCoarseDM(dm, adapt));
462:     break;
463:   case DM_ADAPT_COARSEN:
464:   case DM_ADAPT_COARSEN_LAST:
465:     PetscCall(DMGetFineDM(dm, adapt));
466:     break;
467:   default:
468:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "invalid adaptivity purpose");
469:   }
470:   PetscFunctionReturn(PETSC_SUCCESS);
471: }

473: /*@
474:   DMForestSetAdaptivityPurpose - During the pre-setup phase, set whether the current `DM` is being adapted from its
475:   source (set with `DMForestSetAdaptivityForest()`) for the purpose of refinement (`DM_ADAPT_REFINE`), coarsening
476:   (`DM_ADAPT_COARSEN`), or undefined (`DM_ADAPT_DETERMINE`).

478:   Logically Collective

480:   Input Parameters:
481: + dm      - the forest
482: - purpose - the adaptivity purpose

484:   Level: advanced

486:   Notes:
487:   This only matters for reference counting during `DMDestroy()`. Cyclic references
488:   can be found between `DM`s only if the cyclic reference is due to a fine/coarse relationship
489:   (see `DMSetFineDM()`/`DMSetCoarseDM()`).  If the purpose is not refinement or coarsening, and
490:   the user does not maintain a reference to the post-adaptation forest (i.e., the one created
491:   by `DMForestTemplate()`), this can cause a memory leak.  This method is used by subtypes
492:   of `DMFOREST` when automatically constructing mesh hierarchies.

494: .seealso: `DM`, `DMFOREST`, `DMForestTemplate()`, `DMForestSetAdaptivityForest()`, `DMForestGetAdaptivityForest()`, `DMAdaptFlag`
495: @*/
496: PetscErrorCode DMForestSetAdaptivityPurpose(DM dm, DMAdaptFlag purpose)
497: {
498:   DM_Forest *forest;

500:   PetscFunctionBegin;
501:   forest = (DM_Forest *)dm->data;
502:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the adaptation forest after setup");
503:   if (purpose != forest->adaptPurpose) {
504:     DM adapt;

506:     PetscCall(DMForestGetAdaptivityForest(dm, &adapt));
507:     PetscCall(PetscObjectReference((PetscObject)adapt));
508:     PetscCall(DMForestSetAdaptivityForest(dm, NULL));

510:     forest->adaptPurpose = purpose;

512:     PetscCall(DMForestSetAdaptivityForest(dm, adapt));
513:     PetscCall(DMDestroy(&adapt));
514:   }
515:   PetscFunctionReturn(PETSC_SUCCESS);
516: }

518: /*@
519:   DMForestGetAdaptivityPurpose - Get whether the current `DM` is being adapted from its source (set with
520:   `DMForestSetAdaptivityForest()`) for the purpose of refinement (`DM_ADAPT_REFINE`), coarsening (`DM_ADAPT_COARSEN`),
521:   coarsening only the last level (`DM_ADAPT_COARSEN_LAST`) or undefined (`DM_ADAPT_DETERMINE`).

523:   Not Collective

525:   Input Parameter:
526: . dm - the forest

528:   Output Parameter:
529: . purpose - the adaptivity purpose

531:   Level: advanced

533:   Notes:
534:   This only matters for reference counting: during `DMDestroy()`. Cyclic references
535:   can be found between `DM`s only if the cyclic reference is due to a fine/coarse relationship
536:   (See `DMSetFineDM()`/`DMSetCoarseDM()`).  If the purpose is not refinement or coarsening, and
537:   the user does not maintain a reference to the post-adaptation forest (i.e., the one created
538:   by `DMForestTemplate()`), this can cause a memory leak.  This method is used by subtypes
539:   of `DMFOREST` when automatically constructing mesh hierarchies.

541: .seealso: `DM`, `DMFOREST`, `DMForestTemplate()`, `DMForestSetAdaptivityForest()`, `DMForestGetAdaptivityForest()`, `DMAdaptFlag`
542: @*/
543: PetscErrorCode DMForestGetAdaptivityPurpose(DM dm, DMAdaptFlag *purpose)
544: {
545:   DM_Forest *forest;

547:   PetscFunctionBegin;
548:   forest   = (DM_Forest *)dm->data;
549:   *purpose = forest->adaptPurpose;
550:   PetscFunctionReturn(PETSC_SUCCESS);
551: }

553: /*@
554:   DMForestSetAdjacencyDimension - During the pre-setup phase, set the dimension of interface points that determine
555:   cell adjacency (for the purposes of partitioning and overlap).

557:   Logically Collective

559:   Input Parameters:
560: + dm     - the forest
561: - adjDim - default 0 (i.e., vertices determine adjacency)

563:   Level: intermediate

565: .seealso: `DM`, `DMFOREST`, `DMForestGetAdjacencyDimension()`, `DMForestSetAdjacencyCodimension()`, `DMForestSetPartitionOverlap()`
566: @*/
567: PetscErrorCode DMForestSetAdjacencyDimension(DM dm, PetscInt adjDim)
568: {
569:   PetscInt   dim;
570:   DM_Forest *forest = (DM_Forest *)dm->data;

572:   PetscFunctionBegin;
574:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the adjacency dimension after setup");
575:   PetscCheck(adjDim >= 0, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "adjacency dim cannot be < 0: %" PetscInt_FMT, adjDim);
576:   PetscCall(DMGetDimension(dm, &dim));
577:   PetscCheck(adjDim <= dim, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "adjacency dim cannot be > %" PetscInt_FMT ": %" PetscInt_FMT, dim, adjDim);
578:   forest->adjDim = adjDim;
579:   PetscFunctionReturn(PETSC_SUCCESS);
580: }

582: /*@
583:   DMForestSetAdjacencyCodimension - Like `DMForestSetAdjacencyDimension()`, but specified as a co-dimension (so that,
584:   e.g., adjacency based on facets can be specified by codimension 1 in all cases)

586:   Logically Collective

588:   Input Parameters:
589: + dm       - the forest
590: - adjCodim - default is the dimension of the forest (see `DMGetDimension()`), since this is the codimension of vertices

592:   Level: intermediate

594: .seealso: `DM`, `DMFOREST`, `DMForestGetAdjacencyCodimension()`, `DMForestSetAdjacencyDimension()`
595: @*/
596: PetscErrorCode DMForestSetAdjacencyCodimension(DM dm, PetscInt adjCodim)
597: {
598:   PetscInt dim;

600:   PetscFunctionBegin;
602:   PetscCall(DMGetDimension(dm, &dim));
603:   PetscCall(DMForestSetAdjacencyDimension(dm, dim - adjCodim));
604:   PetscFunctionReturn(PETSC_SUCCESS);
605: }

607: /*@
608:   DMForestGetAdjacencyDimension - Get the dimension of interface points that determine cell adjacency (for the
609:   purposes of partitioning and overlap).

611:   Not Collective

613:   Input Parameter:
614: . dm - the forest

616:   Output Parameter:
617: . adjDim - default 0 (i.e., vertices determine adjacency)

619:   Level: intermediate

621: .seealso: `DM`, `DMFOREST`, `DMForestSetAdjacencyDimension()`, `DMForestGetAdjacencyCodimension()`, `DMForestSetPartitionOverlap()`
622: @*/
623: PetscErrorCode DMForestGetAdjacencyDimension(DM dm, PetscInt *adjDim)
624: {
625:   DM_Forest *forest = (DM_Forest *)dm->data;

627:   PetscFunctionBegin;
629:   PetscAssertPointer(adjDim, 2);
630:   *adjDim = forest->adjDim;
631:   PetscFunctionReturn(PETSC_SUCCESS);
632: }

634: /*@
635:   DMForestGetAdjacencyCodimension - Like `DMForestGetAdjacencyDimension()`, but specified as a co-dimension (so that,
636:   e.g., adjacency based on facets can be specified by codimension 1 in all cases)

638:   Not Collective

640:   Input Parameter:
641: . dm - the forest

643:   Output Parameter:
644: . adjCodim - default isthe dimension of the forest (see `DMGetDimension()`), since this is the codimension of vertices

646:   Level: intermediate

648: .seealso: `DM`, `DMFOREST`, `DMForestSetAdjacencyCodimension()`, `DMForestGetAdjacencyDimension()`
649: @*/
650: PetscErrorCode DMForestGetAdjacencyCodimension(DM dm, PetscInt *adjCodim)
651: {
652:   DM_Forest *forest = (DM_Forest *)dm->data;
653:   PetscInt   dim;

655:   PetscFunctionBegin;
657:   PetscAssertPointer(adjCodim, 2);
658:   PetscCall(DMGetDimension(dm, &dim));
659:   *adjCodim = dim - forest->adjDim;
660:   PetscFunctionReturn(PETSC_SUCCESS);
661: }

663: /*@
664:   DMForestSetPartitionOverlap - During the pre-setup phase, set the amount of cell-overlap present in parallel
665:   partitions of a forest, with values > 0 indicating subdomains that are expanded by that many iterations of adding
666:   adjacent cells

668:   Logically Collective

670:   Input Parameters:
671: + dm      - the forest
672: - overlap - default 0

674:   Level: intermediate

676: .seealso: `DM`, `DMFOREST`, `DMForestGetPartitionOverlap()`, `DMForestSetAdjacencyDimension()`, `DMForestSetAdjacencyCodimension()`
677: @*/
678: PetscErrorCode DMForestSetPartitionOverlap(DM dm, PetscInt overlap)
679: {
680:   DM_Forest *forest = (DM_Forest *)dm->data;

682:   PetscFunctionBegin;
684:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the overlap after setup");
685:   PetscCheck(overlap >= 0, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "overlap cannot be < 0: %" PetscInt_FMT, overlap);
686:   forest->overlap = overlap;
687:   PetscFunctionReturn(PETSC_SUCCESS);
688: }

690: /*@
691:   DMForestGetPartitionOverlap - Get the amount of cell-overlap present in parallel partitions of a forest, with values
692:   > 0 indicating subdomains that are expanded by that many iterations of adding adjacent cells

694:   Not Collective

696:   Input Parameter:
697: . dm - the forest

699:   Output Parameter:
700: . overlap - default 0

702:   Level: intermediate

704: .seealso: `DM`, `DMFOREST`, `DMForestSetAdjacencyDimension()`, `DMForestSetAdjacencyCodimension()`
705: @*/
706: PetscErrorCode DMForestGetPartitionOverlap(DM dm, PetscInt *overlap)
707: {
708:   DM_Forest *forest = (DM_Forest *)dm->data;

710:   PetscFunctionBegin;
712:   PetscAssertPointer(overlap, 2);
713:   *overlap = forest->overlap;
714:   PetscFunctionReturn(PETSC_SUCCESS);
715: }

717: /*@
718:   DMForestSetMinimumRefinement - During the pre-setup phase, set the minimum level of refinement (relative to the base
719:   `DM`, see `DMForestGetBaseDM()`) allowed in the forest.  If the forest is being created by coarsening a previous forest
720:   (see `DMForestGetAdaptivityForest()`) this limits the amount of coarsening.

722:   Logically Collective

724:   Input Parameters:
725: + dm            - the forest
726: - minRefinement - default `PETSC_DEFAULT` (interpreted by the subtype of `DMFOREST`)

728:   Level: intermediate

730: .seealso: `DM`, `DMFOREST`, `DMForestGetMinimumRefinement()`, `DMForestSetMaximumRefinement()`, `DMForestSetInitialRefinement()`, `DMForestGetBaseDM()`, `DMForestGetAdaptivityForest()`
731: @*/
732: PetscErrorCode DMForestSetMinimumRefinement(DM dm, PetscInt minRefinement)
733: {
734:   DM_Forest *forest = (DM_Forest *)dm->data;

736:   PetscFunctionBegin;
738:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the minimum refinement after setup");
739:   forest->minRefinement = minRefinement;
740:   PetscFunctionReturn(PETSC_SUCCESS);
741: }

743: /*@
744:   DMForestGetMinimumRefinement - Get the minimum level of refinement (relative to the base `DM`, see
745:   `DMForestGetBaseDM()`) allowed in the forest.  If the forest is being created by coarsening a previous forest (see
746:   `DMForestGetAdaptivityForest()`), this limits the amount of coarsening.

748:   Not Collective

750:   Input Parameter:
751: . dm - the forest

753:   Output Parameter:
754: . minRefinement - default `PETSC_DEFAULT` (interpreted by the subtype of `DMFOREST`)

756:   Level: intermediate

758: .seealso: `DM`, `DMFOREST`, `DMForestSetMinimumRefinement()`, `DMForestGetMaximumRefinement()`, `DMForestGetInitialRefinement()`, `DMForestGetBaseDM()`, `DMForestGetAdaptivityForest()`
759: @*/
760: PetscErrorCode DMForestGetMinimumRefinement(DM dm, PetscInt *minRefinement)
761: {
762:   DM_Forest *forest = (DM_Forest *)dm->data;

764:   PetscFunctionBegin;
766:   PetscAssertPointer(minRefinement, 2);
767:   *minRefinement = forest->minRefinement;
768:   PetscFunctionReturn(PETSC_SUCCESS);
769: }

771: /*@
772:   DMForestSetInitialRefinement - During the pre-setup phase, set the initial level of refinement (relative to the base
773:   `DM`, see `DMForestGetBaseDM()`) allowed in the forest.

775:   Logically Collective

777:   Input Parameters:
778: + dm             - the forest
779: - initRefinement - default `PETSC_DEFAULT` (interpreted by the subtype of `DMFOREST`)

781:   Level: intermediate

783: .seealso: `DM`, `DMFOREST`, `DMForestSetMinimumRefinement()`, `DMForestSetMaximumRefinement()`, `DMForestGetBaseDM()`
784: @*/
785: PetscErrorCode DMForestSetInitialRefinement(DM dm, PetscInt initRefinement)
786: {
787:   DM_Forest *forest = (DM_Forest *)dm->data;

789:   PetscFunctionBegin;
791:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the initial refinement after setup");
792:   forest->initRefinement = initRefinement;
793:   PetscFunctionReturn(PETSC_SUCCESS);
794: }

796: /*@
797:   DMForestGetInitialRefinement - Get the initial level of refinement (relative to the base `DM`, see
798:   `DMForestGetBaseDM()`) allowed in the forest.

800:   Not Collective

802:   Input Parameter:
803: . dm - the forest

805:   Output Parameter:
806: . initRefinement - default `PETSC_DEFAULT` (interpreted by the subtype of `DMFOREST`)

808:   Level: intermediate

810: .seealso: `DM`, `DMFOREST`, `DMForestSetMinimumRefinement()`, `DMForestSetMaximumRefinement()`, `DMForestGetBaseDM()`
811: @*/
812: PetscErrorCode DMForestGetInitialRefinement(DM dm, PetscInt *initRefinement)
813: {
814:   DM_Forest *forest = (DM_Forest *)dm->data;

816:   PetscFunctionBegin;
818:   PetscAssertPointer(initRefinement, 2);
819:   *initRefinement = forest->initRefinement;
820:   PetscFunctionReturn(PETSC_SUCCESS);
821: }

823: /*@
824:   DMForestSetMaximumRefinement - During the pre-setup phase, set the maximum level of refinement (relative to the base
825:   `DM`, see `DMForestGetBaseDM()`) allowed in the forest.  If the forest is being created by refining a previous forest
826:   (see `DMForestGetAdaptivityForest()`), this limits the amount of refinement.

828:   Logically Collective

830:   Input Parameters:
831: + dm            - the forest
832: - maxRefinement - default `PETSC_DEFAULT` (interpreted by the subtype of `DMFOREST`)

834:   Level: intermediate

836: .seealso: `DM`, `DMFOREST`, `DMForestGetMinimumRefinement()`, `DMForestSetInitialRefinement()`, `DMForestGetBaseDM()`, `DMForestGetAdaptivityDM()`
837: @*/
838: PetscErrorCode DMForestSetMaximumRefinement(DM dm, PetscInt maxRefinement)
839: {
840:   DM_Forest *forest = (DM_Forest *)dm->data;

842:   PetscFunctionBegin;
844:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the maximum refinement after setup");
845:   forest->maxRefinement = maxRefinement;
846:   PetscFunctionReturn(PETSC_SUCCESS);
847: }

849: /*@
850:   DMForestGetMaximumRefinement - Get the maximum level of refinement (relative to the base `DM`, see
851:   `DMForestGetBaseDM()`) allowed in the forest.  If the forest is being created by refining a previous forest (see
852:   `DMForestGetAdaptivityForest()`), this limits the amount of refinement.

854:   Not Collective

856:   Input Parameter:
857: . dm - the forest

859:   Output Parameter:
860: . maxRefinement - default `PETSC_DEFAULT` (interpreted by the subtype of `DMFOREST`)

862:   Level: intermediate

864: .seealso: `DM`, `DMFOREST`, `DMForestSetMaximumRefinement()`, `DMForestGetMinimumRefinement()`, `DMForestGetInitialRefinement()`, `DMForestGetBaseDM()`, `DMForestGetAdaptivityForest()`
865: @*/
866: PetscErrorCode DMForestGetMaximumRefinement(DM dm, PetscInt *maxRefinement)
867: {
868:   DM_Forest *forest = (DM_Forest *)dm->data;

870:   PetscFunctionBegin;
872:   PetscAssertPointer(maxRefinement, 2);
873:   *maxRefinement = forest->maxRefinement;
874:   PetscFunctionReturn(PETSC_SUCCESS);
875: }

877: /*@
878:   DMForestSetAdaptivityStrategy - During the pre-setup phase, set the strategy for combining adaptivity labels from multiple processes.

880:   Logically Collective

882:   Input Parameters:
883: + dm            - the forest
884: - adaptStrategy - default `DMFORESTADAPTALL`

886:   Level: advanced

888:   Notes:
889:   Subtypes of `DMFOREST` may define their own strategies.  Two default strategies are `DMFORESTADAPTALL`, which indicates that all processes must agree
890:   for a refinement/coarsening flag to be valid, and `DMFORESTADAPTANY`, which indicates that only one process needs to
891:   specify refinement/coarsening.

893: .seealso: `DM`, `DMFOREST`, `DMForestGetAdaptivityStrategy()`, `DMFORESTADAPTALL`, `DMFORESTADAPTANY`
894: @*/
895: PetscErrorCode DMForestSetAdaptivityStrategy(DM dm, DMForestAdaptivityStrategy adaptStrategy)
896: {
897:   DM_Forest *forest = (DM_Forest *)dm->data;

899:   PetscFunctionBegin;
901:   PetscCall(PetscFree(forest->adaptStrategy));
902:   PetscCall(PetscStrallocpy((const char *)adaptStrategy, (char **)&forest->adaptStrategy));
903:   PetscFunctionReturn(PETSC_SUCCESS);
904: }

906: /*@
907:   DMForestGetAdaptivityStrategy - Get the strategy for combining adaptivity labels from multiple processes.

909:   Not Collective

911:   Input Parameter:
912: . dm - the forest

914:   Output Parameter:
915: . adaptStrategy - the adaptivity strategy (default `DMFORESTADAPTALL`)

917:   Level: advanced

919:   Note:
920:   Subtypes
921:   of `DMFOREST` may define their own strategies.  Two default strategies are `DMFORESTADAPTALL`, which indicates that all
922:   processes must agree for a refinement/coarsening flag to be valid, and `DMFORESTADAPTANY`, which indicates that only
923:   one process needs to specify refinement/coarsening.

925: .seealso: `DM`, `DMFOREST`, `DMFORESTADAPTALL`, `DMFORESTADAPTANY`, `DMForestSetAdaptivityStrategy()`
926: @*/
927: PetscErrorCode DMForestGetAdaptivityStrategy(DM dm, DMForestAdaptivityStrategy *adaptStrategy)
928: {
929:   DM_Forest *forest = (DM_Forest *)dm->data;

931:   PetscFunctionBegin;
933:   PetscAssertPointer(adaptStrategy, 2);
934:   *adaptStrategy = forest->adaptStrategy;
935:   PetscFunctionReturn(PETSC_SUCCESS);
936: }

938: /*@
939:   DMForestGetAdaptivitySuccess - Return whether the requested adaptation (refinement, coarsening, repartitioning,
940:   etc.) was successful.

942:   Collective

944:   Input Parameter:
945: . dm - the post-adaptation forest

947:   Output Parameter:
948: . success - `PETSC_TRUE` if the post-adaptation forest is different from the pre-adaptation forest.

950:   Level: intermediate

952:   Notes:
953:   `PETSC_FALSE` indicates that the post-adaptation forest is the same as the pre-adaptation
954:   forest.  A requested adaptation may have been unsuccessful if, for example, the requested refinement would have
955:   exceeded the maximum refinement level.

957: .seealso: `DM`, `DMFOREST`
958: @*/
959: PetscErrorCode DMForestGetAdaptivitySuccess(DM dm, PetscBool *success)
960: {
961:   DM_Forest *forest;

963:   PetscFunctionBegin;
965:   PetscCheck(dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "DMSetUp() has not been called yet.");
966:   forest = (DM_Forest *)dm->data;
967:   PetscCall(forest->getadaptivitysuccess(dm, success));
968:   PetscFunctionReturn(PETSC_SUCCESS);
969: }

971: /*@
972:   DMForestSetComputeAdaptivitySF - During the pre-setup phase, set whether transfer `PetscSF`s should be computed
973:   relating the cells of the pre-adaptation forest to the post-adaptiation forest.

975:   Logically Collective

977:   Input Parameters:
978: + dm        - the post-adaptation forest
979: - computeSF - default `PETSC_TRUE`

981:   Level: advanced

983:   Note:
984:   After `DMSetUp()` is called, the transfer `PetscSF`s can be accessed with `DMForestGetAdaptivitySF()`.

986: .seealso: `DM`, `DMFOREST`, `DMForestGetComputeAdaptivitySF()`, `DMForestGetAdaptivitySF()`
987: @*/
988: PetscErrorCode DMForestSetComputeAdaptivitySF(DM dm, PetscBool computeSF)
989: {
990:   DM_Forest *forest;

992:   PetscFunctionBegin;
994:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot compute adaptivity PetscSFs after setup is called");
995:   forest                 = (DM_Forest *)dm->data;
996:   forest->computeAdaptSF = computeSF;
997:   PetscFunctionReturn(PETSC_SUCCESS);
998: }

1000: PetscErrorCode DMForestTransferVec(DM dmIn, Vec vecIn, DM dmOut, Vec vecOut, PetscBool useBCs, PetscReal time)
1001: {
1002:   DM_Forest *forest;

1004:   PetscFunctionBegin;
1009:   forest = (DM_Forest *)dmIn->data;
1010:   PetscCheck(forest->transfervec, PetscObjectComm((PetscObject)dmIn), PETSC_ERR_SUP, "DMForestTransferVec() not implemented");
1011:   PetscCall(forest->transfervec(dmIn, vecIn, dmOut, vecOut, useBCs, time));
1012:   PetscFunctionReturn(PETSC_SUCCESS);
1013: }

1015: PetscErrorCode DMForestTransferVecFromBase(DM dm, Vec vecIn, Vec vecOut)
1016: {
1017:   DM_Forest *forest;

1019:   PetscFunctionBegin;
1023:   forest = (DM_Forest *)dm->data;
1024:   PetscCheck(forest->transfervecfrombase, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "DMForestTransferVecFromBase() not implemented");
1025:   PetscCall(forest->transfervecfrombase(dm, vecIn, vecOut));
1026:   PetscFunctionReturn(PETSC_SUCCESS);
1027: }

1029: /*@
1030:   DMForestGetComputeAdaptivitySF - Get whether transfer `PetscSF`s should be computed relating the cells of the
1031:   pre-adaptation forest to the post-adaptiation forest.  After `DMSetUp()` is called, these transfer PetscSFs can be
1032:   accessed with `DMForestGetAdaptivitySF()`.

1034:   Not Collective

1036:   Input Parameter:
1037: . dm - the post-adaptation forest

1039:   Output Parameter:
1040: . computeSF - default `PETSC_TRUE`

1042:   Level: advanced

1044: .seealso: `DM`, `DMFOREST`, `DMForestSetComputeAdaptivitySF()`, `DMForestGetAdaptivitySF()`
1045: @*/
1046: PetscErrorCode DMForestGetComputeAdaptivitySF(DM dm, PetscBool *computeSF)
1047: {
1048:   DM_Forest *forest;

1050:   PetscFunctionBegin;
1052:   forest     = (DM_Forest *)dm->data;
1053:   *computeSF = forest->computeAdaptSF;
1054:   PetscFunctionReturn(PETSC_SUCCESS);
1055: }

1057: /*@
1058:   DMForestGetAdaptivitySF - Get `PetscSF`s that relate the pre-adaptation forest to the
1059:   post-adaptation forest.

1061:   Not Collective

1063:   Input Parameter:
1064: . dm - the post-adaptation forest

1066:   Output Parameters:
1067: + preCoarseToFine - pre-adaptation coarse cells to post-adaptation fine cells: BCast goes from pre- to post-
1068: - coarseToPreFine - post-adaptation coarse cells to pre-adaptation fine cells: BCast goes from post- to pre-

1070:   Level: advanced

1072:   Notes:
1073:   Adaptation can be any combination of refinement, coarsening, repartition, and change of
1074:   overlap, so there may be some cells of the pre-adaptation that are parents of post-adaptation
1075:   cells, and vice versa.  Therefore there are two `PetscSF`s: one that relates pre-adaptation
1076:   coarse cells to post-adaptation fine cells, and one that relates pre-adaptation fine cells to
1077:   post-adaptation coarse cells.

1079: .seealso: `DM`, `DMFOREST`, `DMForestGetComputeAdaptivitySF()`, `DMForestSetComputeAdaptivitySF()`
1080: @*/
1081: PetscErrorCode DMForestGetAdaptivitySF(DM dm, PetscSF *preCoarseToFine, PetscSF *coarseToPreFine)
1082: {
1083:   DM_Forest *forest;

1085:   PetscFunctionBegin;
1087:   PetscCall(DMSetUp(dm));
1088:   forest = (DM_Forest *)dm->data;
1089:   if (preCoarseToFine) *preCoarseToFine = forest->preCoarseToFine;
1090:   if (coarseToPreFine) *coarseToPreFine = forest->coarseToPreFine;
1091:   PetscFunctionReturn(PETSC_SUCCESS);
1092: }

1094: /*@
1095:   DMForestSetGradeFactor - During the pre-setup phase, set the desired amount of grading in the
1096:   mesh, e.g. give 2 to indicate that the diameter of neighboring cells should differ by at most
1097:   a factor of 2.

1099:   Logically Collective

1101:   Input Parameters:
1102: + dm    - the forest
1103: - grade - the grading factor

1105:   Level: advanced

1107:   Note:
1108:   Subtypes of `DMFOREST` may only support one particular choice of grading factor.

1110: .seealso: `DM`, `DMFOREST`, `DMForestGetGradeFactor()`
1111: @*/
1112: PetscErrorCode DMForestSetGradeFactor(DM dm, PetscInt grade)
1113: {
1114:   DM_Forest *forest = (DM_Forest *)dm->data;

1116:   PetscFunctionBegin;
1118:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the grade factor after setup");
1119:   forest->gradeFactor = grade;
1120:   PetscFunctionReturn(PETSC_SUCCESS);
1121: }

1123: /*@
1124:   DMForestGetGradeFactor - Get the desired amount of grading in the mesh, e.g. give 2 to indicate that the diameter of
1125:   neighboring cells should differ by at most a factor of 2.  Subtypes of `DMFOREST` may only support one particular
1126:   choice of grading factor.

1128:   Not Collective

1130:   Input Parameter:
1131: . dm - the forest

1133:   Output Parameter:
1134: . grade - the grading factor

1136:   Level: advanced

1138: .seealso: `DM`, `DMFOREST`, `DMForestSetGradeFactor()`
1139: @*/
1140: PetscErrorCode DMForestGetGradeFactor(DM dm, PetscInt *grade)
1141: {
1142:   DM_Forest *forest = (DM_Forest *)dm->data;

1144:   PetscFunctionBegin;
1146:   PetscAssertPointer(grade, 2);
1147:   *grade = forest->gradeFactor;
1148:   PetscFunctionReturn(PETSC_SUCCESS);
1149: }

1151: /*@
1152:   DMForestSetCellWeightFactor - During the pre-setup phase, set the factor by which the level of refinement changes
1153:   the cell weight (see `DMForestSetCellWeights()`) when calculating partitions.

1155:   Logically Collective

1157:   Input Parameters:
1158: + dm            - the forest
1159: - weightsFactor - default 1.

1161:   Level: advanced

1163:   Note:
1164:   The final weight of a cell will be (cellWeight) * (weightFactor^refinementLevel).  A factor
1165:   of 1 indicates that the weight of a cell does not depend on its level; a factor of 2, for
1166:   example, might be appropriate for sub-cycling time-stepping methods, when the computation
1167:   associated with a cell is multiplied by a factor of 2 for each additional level of
1168:   refinement.

1170: .seealso: `DM`, `DMFOREST`, `DMForestGetCellWeightFactor()`, `DMForestSetCellWeights()`
1171: @*/
1172: PetscErrorCode DMForestSetCellWeightFactor(DM dm, PetscReal weightsFactor)
1173: {
1174:   DM_Forest *forest = (DM_Forest *)dm->data;

1176:   PetscFunctionBegin;
1178:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the weights factor after setup");
1179:   forest->weightsFactor = weightsFactor;
1180:   PetscFunctionReturn(PETSC_SUCCESS);
1181: }

1183: /*@
1184:   DMForestGetCellWeightFactor - Get the factor by which the level of refinement changes the cell weight (see
1185:   `DMForestSetCellWeights()`) when calculating partitions.

1187:   Not Collective

1189:   Input Parameter:
1190: . dm - the forest

1192:   Output Parameter:
1193: . weightsFactor - default 1.

1195:   Level: advanced

1197:   Note:
1198:   The final weight of a cell will be (cellWeight) * (weightFactor^refinementLevel).  A factor
1199:   of 1 indicates that the weight of a cell does not depend on its level; a factor of 2, for
1200:   example, might be appropriate for sub-cycling time-stepping methods, when the computation
1201:   associated with a cell is multiplied by a factor of 2 for each additional level of
1202:   refinement.

1204: .seealso: `DM`, `DMFOREST`, `DMForestSetCellWeightFactor()`, `DMForestSetCellWeights()`
1205: @*/
1206: PetscErrorCode DMForestGetCellWeightFactor(DM dm, PetscReal *weightsFactor)
1207: {
1208:   DM_Forest *forest = (DM_Forest *)dm->data;

1210:   PetscFunctionBegin;
1212:   PetscAssertPointer(weightsFactor, 2);
1213:   *weightsFactor = forest->weightsFactor;
1214:   PetscFunctionReturn(PETSC_SUCCESS);
1215: }

1217: /*@
1218:   DMForestGetCellChart - After the setup phase, get the local half-open interval of the chart of cells on this process

1220:   Not Collective

1222:   Input Parameter:
1223: . dm - the forest

1225:   Output Parameters:
1226: + cStart - the first cell on this process
1227: - cEnd   - one after the final cell on this process

1229:   Level: intermediate

1231: .seealso: `DM`, `DMFOREST`, `DMForestGetCellSF()`
1232: @*/
1233: PetscErrorCode DMForestGetCellChart(DM dm, PetscInt *cStart, PetscInt *cEnd)
1234: {
1235:   DM_Forest *forest = (DM_Forest *)dm->data;

1237:   PetscFunctionBegin;
1239:   PetscAssertPointer(cStart, 2);
1240:   PetscAssertPointer(cEnd, 3);
1241:   if (((forest->cStart == PETSC_DETERMINE) || (forest->cEnd == PETSC_DETERMINE)) && forest->createcellchart) PetscCall(forest->createcellchart(dm, &forest->cStart, &forest->cEnd));
1242:   *cStart = forest->cStart;
1243:   *cEnd   = forest->cEnd;
1244:   PetscFunctionReturn(PETSC_SUCCESS);
1245: }

1247: /*@
1248:   DMForestGetCellSF - After the setup phase, get the `PetscSF` for overlapping cells between processes

1250:   Not Collective

1252:   Input Parameter:
1253: . dm - the forest

1255:   Output Parameter:
1256: . cellSF - the `PetscSF`

1258:   Level: intermediate

1260: .seealso: `DM`, `DMFOREST`, `DMForestGetCellChart()`
1261: @*/
1262: PetscErrorCode DMForestGetCellSF(DM dm, PetscSF *cellSF)
1263: {
1264:   DM_Forest *forest = (DM_Forest *)dm->data;

1266:   PetscFunctionBegin;
1268:   PetscAssertPointer(cellSF, 2);
1269:   if ((!forest->cellSF) && forest->createcellsf) PetscCall(forest->createcellsf(dm, &forest->cellSF));
1270:   *cellSF = forest->cellSF;
1271:   PetscFunctionReturn(PETSC_SUCCESS);
1272: }

1274: /*@
1275:   DMForestSetAdaptivityLabel - During the pre-setup phase, set the label of the pre-adaptation forest (see
1276:   `DMForestGetAdaptivityForest()`) that holds the adaptation flags (refinement, coarsening, or some combination).

1278:   Logically Collective

1280:   Input Parameters:
1281: + dm         - the forest
1282: - adaptLabel - the label in the pre-adaptation forest

1284:   Level: intermediate

1286:   Note:
1287:   The interpretation of the label values is up to the subtype of `DMFOREST`, but
1288:   `DM_ADAPT_DETERMINE`, `DM_ADAPT_KEEP`, `DM_ADAPT_REFINE`, and `DM_ADAPT_COARSEN` have been
1289:   reserved as choices that should be accepted by all subtypes.

1291: .seealso: `DM`, `DMFOREST`, `DMForestGetAdaptivityLabel()`
1292: @*/
1293: PetscErrorCode DMForestSetAdaptivityLabel(DM dm, DMLabel adaptLabel)
1294: {
1295:   DM_Forest *forest = (DM_Forest *)dm->data;

1297:   PetscFunctionBegin;
1300:   PetscCall(PetscObjectReference((PetscObject)adaptLabel));
1301:   PetscCall(DMLabelDestroy(&forest->adaptLabel));
1302:   forest->adaptLabel = adaptLabel;
1303:   PetscFunctionReturn(PETSC_SUCCESS);
1304: }

1306: /*@
1307:   DMForestGetAdaptivityLabel - Get the label of the pre-adaptation forest (see `DMForestGetAdaptivityForest()`) that
1308:   holds the adaptation flags (refinement, coarsening, or some combination).

1310:   Not Collective

1312:   Input Parameter:
1313: . dm - the forest

1315:   Output Parameter:
1316: . adaptLabel - the name of the label in the pre-adaptation forest

1318:   Level: intermediate

1320:   Note:
1321:   The interpretation of the label values is up to the subtype of `DMFOREST`, but
1322:   `DM_ADAPT_DETERMINE`, `DM_ADAPT_KEEP`, `DM_ADAPT_REFINE`, and `DM_ADAPT_COARSEN` have been
1323:   reserved as choices that should be accepted by all subtypes.

1325: .seealso: `DM`, `DMFOREST`, `DMForestSetAdaptivityLabel()`
1326: @*/
1327: PetscErrorCode DMForestGetAdaptivityLabel(DM dm, DMLabel *adaptLabel)
1328: {
1329:   DM_Forest *forest = (DM_Forest *)dm->data;

1331:   PetscFunctionBegin;
1333:   *adaptLabel = forest->adaptLabel;
1334:   PetscFunctionReturn(PETSC_SUCCESS);
1335: }

1337: /*@
1338:   DMForestSetCellWeights - Set the weights assigned to each of the cells (see `DMForestGetCellChart()`) of the current
1339:   process: weights are used to determine parallel partitioning.

1341:   Logically Collective

1343:   Input Parameters:
1344: + dm       - the forest
1345: . weights  - the array of weights (see `DMForestSetWeightCapacity()`) for all cells, or `NULL` to indicate each cell has weight 1.
1346: - copyMode - how weights should reference weights

1348:   Level: advanced

1350: .seealso: `DM`, `DMFOREST`, `DMForestGetCellWeights()`, `DMForestSetWeightCapacity()`
1351: @*/
1352: PetscErrorCode DMForestSetCellWeights(DM dm, PetscReal weights[], PetscCopyMode copyMode)
1353: {
1354:   DM_Forest *forest = (DM_Forest *)dm->data;
1355:   PetscInt   cStart, cEnd;

1357:   PetscFunctionBegin;
1359:   PetscCall(DMForestGetCellChart(dm, &cStart, &cEnd));
1360:   PetscCheck(cEnd >= cStart, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "cell chart [%" PetscInt_FMT ",%" PetscInt_FMT ") is not valid", cStart, cEnd);
1361:   if (copyMode == PETSC_COPY_VALUES) {
1362:     if (forest->cellWeightsCopyMode != PETSC_OWN_POINTER || forest->cellWeights == weights) PetscCall(PetscMalloc1(cEnd - cStart, &forest->cellWeights));
1363:     PetscCall(PetscArraycpy(forest->cellWeights, weights, cEnd - cStart));
1364:     forest->cellWeightsCopyMode = PETSC_OWN_POINTER;
1365:     PetscFunctionReturn(PETSC_SUCCESS);
1366:   }
1367:   if (forest->cellWeightsCopyMode == PETSC_OWN_POINTER) PetscCall(PetscFree(forest->cellWeights));
1368:   forest->cellWeights         = weights;
1369:   forest->cellWeightsCopyMode = copyMode;
1370:   PetscFunctionReturn(PETSC_SUCCESS);
1371: }

1373: /*@
1374:   DMForestGetCellWeights - Get the weights assigned to each of the cells (see `DMForestGetCellChart()`) of the current
1375:   process: weights are used to determine parallel partitioning.

1377:   Not Collective

1379:   Input Parameter:
1380: . dm - the forest

1382:   Output Parameter:
1383: . weights - the array of weights for all cells, or `NULL` to indicate each cell has weight 1.

1385:   Level: advanced

1387: .seealso: `DM`, `DMFOREST`, `DMForestSetCellWeights()`, `DMForestSetWeightCapacity()`
1388: @*/
1389: PetscErrorCode DMForestGetCellWeights(DM dm, PetscReal **weights)
1390: {
1391:   DM_Forest *forest = (DM_Forest *)dm->data;

1393:   PetscFunctionBegin;
1395:   PetscAssertPointer(weights, 2);
1396:   *weights = forest->cellWeights;
1397:   PetscFunctionReturn(PETSC_SUCCESS);
1398: }

1400: /*@
1401:   DMForestSetWeightCapacity - During the pre-setup phase, set the capacity of the current process when repartitioning
1402:   a pre-adaptation forest (see `DMForestGetAdaptivityForest()`).

1404:   Logically Collective

1406:   Input Parameters:
1407: + dm       - the forest
1408: - capacity - this process's capacity

1410:   Level: advanced

1412:   Note:
1413:   After partitioning, the ratio of the weight of each process's cells to the process's capacity
1414:   will be roughly equal for all processes.  A capacity of 0 indicates that the current process
1415:   should not have any cells after repartitioning.

1417: .seealso: `DM`, `DMFOREST`, `DMForestGetWeightCapacity()`, `DMForestSetCellWeights()`, `DMForestSetCellWeightFactor()`
1418: @*/
1419: PetscErrorCode DMForestSetWeightCapacity(DM dm, PetscReal capacity)
1420: {
1421:   DM_Forest *forest = (DM_Forest *)dm->data;

1423:   PetscFunctionBegin;
1425:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the weight capacity after setup");
1426:   PetscCheck(capacity >= 0., PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "Cannot have negative weight capacity; %g", (double)capacity);
1427:   forest->weightCapacity = capacity;
1428:   PetscFunctionReturn(PETSC_SUCCESS);
1429: }

1431: /*@
1432:   DMForestGetWeightCapacity - Set the capacity of the current process when repartitioning a pre-adaptation forest (see
1433:   `DMForestGetAdaptivityForest()`).

1435:   Not Collective

1437:   Input Parameter:
1438: . dm - the forest

1440:   Output Parameter:
1441: . capacity - this process's capacity

1443:   Level: advanced

1445:   Note:
1446:   After partitioning, the ratio of the weight of each process's cells to the process's capacity
1447:   will be roughly equal for all processes.  A capacity of 0 indicates that the current process
1448:   should not have any cells after repartitioning.

1450: .seealso: `DM`, `DMFOREST`, `DMForestSetWeightCapacity()`, `DMForestSetCellWeights()`, `DMForestSetCellWeightFactor()`
1451: @*/
1452: PetscErrorCode DMForestGetWeightCapacity(DM dm, PetscReal *capacity)
1453: {
1454:   DM_Forest *forest = (DM_Forest *)dm->data;

1456:   PetscFunctionBegin;
1458:   PetscAssertPointer(capacity, 2);
1459:   *capacity = forest->weightCapacity;
1460:   PetscFunctionReturn(PETSC_SUCCESS);
1461: }

1463: PETSC_EXTERN PetscErrorCode DMSetFromOptions_Forest(DM dm, PetscOptionItems *PetscOptionsObject)
1464: {
1465:   PetscBool                  flg, flg1, flg2, flg3, flg4;
1466:   DMForestTopology           oldTopo;
1467:   char                       stringBuffer[256];
1468:   PetscViewer                viewer;
1469:   PetscViewerFormat          format;
1470:   PetscInt                   adjDim, adjCodim, overlap, minRefinement, initRefinement, maxRefinement, grade;
1471:   PetscReal                  weightsFactor;
1472:   DMForestAdaptivityStrategy adaptStrategy;

1474:   PetscFunctionBegin;
1475:   PetscCall(DMForestGetTopology(dm, &oldTopo));
1476:   PetscOptionsHeadBegin(PetscOptionsObject, "DMForest Options");
1477:   PetscCall(PetscOptionsString("-dm_forest_topology", "the topology of the forest's base mesh", "DMForestSetTopology", oldTopo, stringBuffer, sizeof(stringBuffer), &flg1));
1478:   PetscCall(PetscOptionsViewer("-dm_forest_base_dm", "load the base DM from a viewer specification", "DMForestSetBaseDM", &viewer, &format, &flg2));
1479:   PetscCall(PetscOptionsViewer("-dm_forest_coarse_forest", "load the coarse forest from a viewer specification", "DMForestSetCoarseForest", &viewer, &format, &flg3));
1480:   PetscCall(PetscOptionsViewer("-dm_forest_fine_forest", "load the fine forest from a viewer specification", "DMForestSetFineForest", &viewer, &format, &flg4));
1481:   PetscCheck((PetscInt)flg1 + (PetscInt)flg2 + (PetscInt)flg3 + (PetscInt)flg4 <= 1, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_INCOMP, "Specify only one of -dm_forest_{topology,base_dm,coarse_forest,fine_forest}");
1482:   if (flg1) {
1483:     PetscCall(DMForestSetTopology(dm, (DMForestTopology)stringBuffer));
1484:     PetscCall(DMForestSetBaseDM(dm, NULL));
1485:     PetscCall(DMForestSetAdaptivityForest(dm, NULL));
1486:   }
1487:   if (flg2) {
1488:     DM base;

1490:     PetscCall(DMCreate(PetscObjectComm((PetscObject)dm), &base));
1491:     PetscCall(PetscViewerPushFormat(viewer, format));
1492:     PetscCall(DMLoad(base, viewer));
1493:     PetscCall(PetscViewerDestroy(&viewer));
1494:     PetscCall(DMForestSetBaseDM(dm, base));
1495:     PetscCall(DMDestroy(&base));
1496:     PetscCall(DMForestSetTopology(dm, NULL));
1497:     PetscCall(DMForestSetAdaptivityForest(dm, NULL));
1498:   }
1499:   if (flg3) {
1500:     DM coarse;

1502:     PetscCall(DMCreate(PetscObjectComm((PetscObject)dm), &coarse));
1503:     PetscCall(PetscViewerPushFormat(viewer, format));
1504:     PetscCall(DMLoad(coarse, viewer));
1505:     PetscCall(PetscViewerDestroy(&viewer));
1506:     PetscCall(DMForestSetAdaptivityForest(dm, coarse));
1507:     PetscCall(DMDestroy(&coarse));
1508:     PetscCall(DMForestSetTopology(dm, NULL));
1509:     PetscCall(DMForestSetBaseDM(dm, NULL));
1510:   }
1511:   if (flg4) {
1512:     DM fine;

1514:     PetscCall(DMCreate(PetscObjectComm((PetscObject)dm), &fine));
1515:     PetscCall(PetscViewerPushFormat(viewer, format));
1516:     PetscCall(DMLoad(fine, viewer));
1517:     PetscCall(PetscViewerDestroy(&viewer));
1518:     PetscCall(DMForestSetAdaptivityForest(dm, fine));
1519:     PetscCall(DMDestroy(&fine));
1520:     PetscCall(DMForestSetTopology(dm, NULL));
1521:     PetscCall(DMForestSetBaseDM(dm, NULL));
1522:   }
1523:   PetscCall(DMForestGetAdjacencyDimension(dm, &adjDim));
1524:   PetscCall(PetscOptionsBoundedInt("-dm_forest_adjacency_dimension", "set the dimension of points that define adjacency in the forest", "DMForestSetAdjacencyDimension", adjDim, &adjDim, &flg, 0));
1525:   if (flg) {
1526:     PetscCall(DMForestSetAdjacencyDimension(dm, adjDim));
1527:   } else {
1528:     PetscCall(DMForestGetAdjacencyCodimension(dm, &adjCodim));
1529:     PetscCall(PetscOptionsBoundedInt("-dm_forest_adjacency_codimension", "set the codimension of points that define adjacency in the forest", "DMForestSetAdjacencyCodimension", adjCodim, &adjCodim, &flg, 1));
1530:     if (flg) PetscCall(DMForestSetAdjacencyCodimension(dm, adjCodim));
1531:   }
1532:   PetscCall(DMForestGetPartitionOverlap(dm, &overlap));
1533:   PetscCall(PetscOptionsBoundedInt("-dm_forest_partition_overlap", "set the degree of partition overlap", "DMForestSetPartitionOverlap", overlap, &overlap, &flg, 0));
1534:   if (flg) PetscCall(DMForestSetPartitionOverlap(dm, overlap));
1535: #if 0
1536:   PetscCall(PetscOptionsBoundedInt("-dm_refine","equivalent to -dm_forest_set_minimum_refinement and -dm_forest_set_initial_refinement with the same value",NULL,minRefinement,&minRefinement,&flg,0));
1537:   if (flg) {
1538:     PetscCall(DMForestSetMinimumRefinement(dm,minRefinement));
1539:     PetscCall(DMForestSetInitialRefinement(dm,minRefinement));
1540:   }
1541:   PetscCall(PetscOptionsBoundedInt("-dm_refine_hierarchy","equivalent to -dm_forest_set_minimum_refinement 0 and -dm_forest_set_initial_refinement",NULL,initRefinement,&initRefinement,&flg,0));
1542:   if (flg) {
1543:     PetscCall(DMForestSetMinimumRefinement(dm,0));
1544:     PetscCall(DMForestSetInitialRefinement(dm,initRefinement));
1545:   }
1546: #endif
1547:   PetscCall(DMForestGetMinimumRefinement(dm, &minRefinement));
1548:   PetscCall(PetscOptionsBoundedInt("-dm_forest_minimum_refinement", "set the minimum level of refinement in the forest", "DMForestSetMinimumRefinement", minRefinement, &minRefinement, &flg, 0));
1549:   if (flg) PetscCall(DMForestSetMinimumRefinement(dm, minRefinement));
1550:   PetscCall(DMForestGetInitialRefinement(dm, &initRefinement));
1551:   PetscCall(PetscOptionsBoundedInt("-dm_forest_initial_refinement", "set the initial level of refinement in the forest", "DMForestSetInitialRefinement", initRefinement, &initRefinement, &flg, 0));
1552:   if (flg) PetscCall(DMForestSetInitialRefinement(dm, initRefinement));
1553:   PetscCall(DMForestGetMaximumRefinement(dm, &maxRefinement));
1554:   PetscCall(PetscOptionsBoundedInt("-dm_forest_maximum_refinement", "set the maximum level of refinement in the forest", "DMForestSetMaximumRefinement", maxRefinement, &maxRefinement, &flg, 0));
1555:   if (flg) PetscCall(DMForestSetMaximumRefinement(dm, maxRefinement));
1556:   PetscCall(DMForestGetAdaptivityStrategy(dm, &adaptStrategy));
1557:   PetscCall(PetscOptionsString("-dm_forest_adaptivity_strategy", "the forest's adaptivity-flag resolution strategy", "DMForestSetAdaptivityStrategy", adaptStrategy, stringBuffer, sizeof(stringBuffer), &flg));
1558:   if (flg) PetscCall(DMForestSetAdaptivityStrategy(dm, (DMForestAdaptivityStrategy)stringBuffer));
1559:   PetscCall(DMForestGetGradeFactor(dm, &grade));
1560:   PetscCall(PetscOptionsBoundedInt("-dm_forest_grade_factor", "grade factor between neighboring cells", "DMForestSetGradeFactor", grade, &grade, &flg, 0));
1561:   if (flg) PetscCall(DMForestSetGradeFactor(dm, grade));
1562:   PetscCall(DMForestGetCellWeightFactor(dm, &weightsFactor));
1563:   PetscCall(PetscOptionsReal("-dm_forest_cell_weight_factor", "multiplying weight factor for cell refinement", "DMForestSetCellWeightFactor", weightsFactor, &weightsFactor, &flg));
1564:   if (flg) PetscCall(DMForestSetCellWeightFactor(dm, weightsFactor));
1565:   PetscOptionsHeadEnd();
1566:   PetscFunctionReturn(PETSC_SUCCESS);
1567: }

1569: static PetscErrorCode DMCreateSubDM_Forest(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1570: {
1571:   PetscFunctionBegin;
1572:   if (subdm) PetscCall(DMClone(dm, subdm));
1573:   PetscCall(DMCreateSectionSubDM(dm, numFields, fields, NULL, NULL, is, subdm));
1574:   PetscFunctionReturn(PETSC_SUCCESS);
1575: }

1577: static PetscErrorCode DMRefine_Forest(DM dm, MPI_Comm comm, DM *dmRefined)
1578: {
1579:   DMLabel refine;
1580:   DM      fineDM;

1582:   PetscFunctionBegin;
1583:   PetscCall(DMGetFineDM(dm, &fineDM));
1584:   if (fineDM) {
1585:     PetscCall(PetscObjectReference((PetscObject)fineDM));
1586:     *dmRefined = fineDM;
1587:     PetscFunctionReturn(PETSC_SUCCESS);
1588:   }
1589:   PetscCall(DMForestTemplate(dm, comm, dmRefined));
1590:   PetscCall(DMGetLabel(dm, "refine", &refine));
1591:   if (!refine) {
1592:     PetscCall(DMLabelCreate(PETSC_COMM_SELF, "refine", &refine));
1593:     PetscCall(DMLabelSetDefaultValue(refine, DM_ADAPT_REFINE));
1594:   } else PetscCall(PetscObjectReference((PetscObject)refine));
1595:   PetscCall(DMForestSetAdaptivityLabel(*dmRefined, refine));
1596:   PetscCall(DMLabelDestroy(&refine));
1597:   PetscFunctionReturn(PETSC_SUCCESS);
1598: }

1600: static PetscErrorCode DMCoarsen_Forest(DM dm, MPI_Comm comm, DM *dmCoarsened)
1601: {
1602:   DMLabel coarsen;
1603:   DM      coarseDM;

1605:   PetscFunctionBegin;
1606:   if (comm != MPI_COMM_NULL) {
1607:     PetscMPIInt mpiComparison;
1608:     MPI_Comm    dmcomm = PetscObjectComm((PetscObject)dm);

1610:     PetscCallMPI(MPI_Comm_compare(comm, dmcomm, &mpiComparison));
1611:     PetscCheck(mpiComparison == MPI_IDENT || mpiComparison == MPI_CONGRUENT, dmcomm, PETSC_ERR_SUP, "No support for different communicators yet");
1612:   }
1613:   PetscCall(DMGetCoarseDM(dm, &coarseDM));
1614:   if (coarseDM) {
1615:     PetscCall(PetscObjectReference((PetscObject)coarseDM));
1616:     *dmCoarsened = coarseDM;
1617:     PetscFunctionReturn(PETSC_SUCCESS);
1618:   }
1619:   PetscCall(DMForestTemplate(dm, comm, dmCoarsened));
1620:   PetscCall(DMForestSetAdaptivityPurpose(*dmCoarsened, DM_ADAPT_COARSEN));
1621:   PetscCall(DMGetLabel(dm, "coarsen", &coarsen));
1622:   if (!coarsen) {
1623:     PetscCall(DMLabelCreate(PETSC_COMM_SELF, "coarsen", &coarsen));
1624:     PetscCall(DMLabelSetDefaultValue(coarsen, DM_ADAPT_COARSEN));
1625:   } else PetscCall(PetscObjectReference((PetscObject)coarsen));
1626:   PetscCall(DMForestSetAdaptivityLabel(*dmCoarsened, coarsen));
1627:   PetscCall(DMLabelDestroy(&coarsen));
1628:   PetscFunctionReturn(PETSC_SUCCESS);
1629: }

1631: PetscErrorCode DMAdaptLabel_Forest(DM dm, PETSC_UNUSED Vec metric, DMLabel label, PETSC_UNUSED DMLabel rgLabel, DM *adaptedDM)
1632: {
1633:   PetscBool success;

1635:   PetscFunctionBegin;
1636:   PetscCall(DMForestTemplate(dm, PetscObjectComm((PetscObject)dm), adaptedDM));
1637:   PetscCall(DMForestSetAdaptivityLabel(*adaptedDM, label));
1638:   PetscCall(DMSetUp(*adaptedDM));
1639:   PetscCall(DMForestGetAdaptivitySuccess(*adaptedDM, &success));
1640:   if (!success) {
1641:     PetscCall(DMDestroy(adaptedDM));
1642:     *adaptedDM = NULL;
1643:   }
1644:   PetscFunctionReturn(PETSC_SUCCESS);
1645: }

1647: static PetscErrorCode DMInitialize_Forest(DM dm)
1648: {
1649:   PetscFunctionBegin;
1650:   PetscCall(PetscMemzero(dm->ops, sizeof(*dm->ops)));

1652:   dm->ops->clone          = DMClone_Forest;
1653:   dm->ops->setfromoptions = DMSetFromOptions_Forest;
1654:   dm->ops->destroy        = DMDestroy_Forest;
1655:   dm->ops->createsubdm    = DMCreateSubDM_Forest;
1656:   dm->ops->refine         = DMRefine_Forest;
1657:   dm->ops->coarsen        = DMCoarsen_Forest;
1658:   PetscFunctionReturn(PETSC_SUCCESS);
1659: }

1661: /*MC

1663:      DMFOREST = "forest" - A DM object that encapsulates a hierarchically refined mesh.  Forests usually have a base `DM`
1664:   (see `DMForestGetBaseDM()`), from which it is refined.  The refinement and partitioning of forests is considered
1665:   immutable after `DMSetUp()` is called.  To adapt a mesh, one should call `DMForestTemplate()` to create a new mesh that
1666:   will default to being identical to it, specify how that mesh should differ, and then calling `DMSetUp()` on the new
1667:   mesh.

1669:   To specify that a mesh should be refined or coarsened from the previous mesh, a label should be defined on the
1670:   previous mesh whose values indicate which cells should be refined (`DM_ADAPT_REFINE`) or coarsened (`DM_ADAPT_COARSEN`)
1671:   and how (subtypes are free to allow additional values for things like anisotropic refinement).  The label should be
1672:   given to the *new* mesh with `DMForestSetAdaptivityLabel()`.

1674:   Level: advanced

1676: .seealso: `DMType`, `DM`, `DMCreate()`, `DMSetType()`, `DMForestGetBaseDM()`, `DMForestSetBaseDM()`, `DMForestTemplate()`, `DMForestSetAdaptivityLabel()`
1677: M*/

1679: PETSC_EXTERN PetscErrorCode DMCreate_Forest(DM dm)
1680: {
1681:   DM_Forest *forest;

1683:   PetscFunctionBegin;
1685:   PetscCall(PetscNew(&forest));
1686:   dm->dim                     = 0;
1687:   dm->data                    = forest;
1688:   forest->refct               = 1;
1689:   forest->data                = NULL;
1690:   forest->topology            = NULL;
1691:   forest->adapt               = NULL;
1692:   forest->base                = NULL;
1693:   forest->adaptPurpose        = DM_ADAPT_DETERMINE;
1694:   forest->adjDim              = PETSC_DEFAULT;
1695:   forest->overlap             = PETSC_DEFAULT;
1696:   forest->minRefinement       = PETSC_DEFAULT;
1697:   forest->maxRefinement       = PETSC_DEFAULT;
1698:   forest->initRefinement      = PETSC_DEFAULT;
1699:   forest->cStart              = PETSC_DETERMINE;
1700:   forest->cEnd                = PETSC_DETERMINE;
1701:   forest->cellSF              = NULL;
1702:   forest->adaptLabel          = NULL;
1703:   forest->gradeFactor         = 2;
1704:   forest->cellWeights         = NULL;
1705:   forest->cellWeightsCopyMode = PETSC_USE_POINTER;
1706:   forest->weightsFactor       = 1.;
1707:   forest->weightCapacity      = 1.;
1708:   PetscCall(DMForestSetAdaptivityStrategy(dm, DMFORESTADAPTALL));
1709:   PetscCall(DMInitialize_Forest(dm));
1710:   PetscFunctionReturn(PETSC_SUCCESS);
1711: }