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.  The new DM reproduces the configuration
102:   of the source, but is not yet setup, so that the user can then define only the ways that the new DM should differ
103:   (by, e.g., refinement or repartitioning).  The source DM is also set as the adaptivity source DM of the new DM (see
104:   DMForestSetAdaptivityForest()).

106:   Collective on dm

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

112:   Output Parameter:
113: . tdm - the new DM object

115:   Level: intermediate

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

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

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

171: static PetscErrorCode DMInitialize_Forest(DM dm);

173: PETSC_EXTERN PetscErrorCode DMClone_Forest(DM dm, DM *newdm)
174: {
175:   DM_Forest  *forest = (DM_Forest *)dm->data;
176:   const char *type;

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

187: static PetscErrorCode DMDestroy_Forest(DM dm)
188: {
189:   DM_Forest *forest = (DM_Forest *)dm->data;

191:   PetscFunctionBegin;
192:   if (--forest->refct > 0) PetscFunctionReturn(PETSC_SUCCESS);
193:   if (forest->destroy) PetscCall((*forest->destroy)(dm));
194:   PetscCall(PetscSFDestroy(&forest->cellSF));
195:   PetscCall(PetscSFDestroy(&forest->preCoarseToFine));
196:   PetscCall(PetscSFDestroy(&forest->coarseToPreFine));
197:   PetscCall(DMLabelDestroy(&forest->adaptLabel));
198:   PetscCall(PetscFree(forest->adaptStrategy));
199:   PetscCall(DMDestroy(&forest->base));
200:   PetscCall(DMDestroy(&forest->adapt));
201:   PetscCall(PetscFree(forest->topology));
202:   PetscCall(PetscFree(forest));
203:   PetscFunctionReturn(PETSC_SUCCESS);
204: }

206: /*@C
207:   DMForestSetTopology - Set the topology of a DMForest during the pre-setup phase.  The topology is a string (e.g.
208:   "cube", "shell") and can be interpreted by subtypes of DMFOREST) to construct the base DM of a forest during
209:   DMSetUp().

211:   Logically collective on dm

213:   Input parameters:
214: + dm - the forest
215: - topology - the topology of the forest

217:   Level: intermediate

219: .seealso: `DMForestGetTopology()`, `DMForestSetBaseDM()`
220: @*/
221: PetscErrorCode DMForestSetTopology(DM dm, DMForestTopology topology)
222: {
223:   DM_Forest *forest = (DM_Forest *)dm->data;

225:   PetscFunctionBegin;
227:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the topology after setup");
228:   PetscCall(PetscFree(forest->topology));
229:   PetscCall(PetscStrallocpy((const char *)topology, (char **)&forest->topology));
230:   PetscFunctionReturn(PETSC_SUCCESS);
231: }

233: /*@C
234:   DMForestGetTopology - Get a string describing the topology of a DMForest.

236:   Not collective

238:   Input parameter:
239: . dm - the forest

241:   Output parameter:
242: . topology - the topology of the forest (e.g., 'cube', 'shell')

244:   Level: intermediate

246: .seealso: `DMForestSetTopology()`
247: @*/
248: PetscErrorCode DMForestGetTopology(DM dm, DMForestTopology *topology)
249: {
250:   DM_Forest *forest = (DM_Forest *)dm->data;

252:   PetscFunctionBegin;
255:   *topology = forest->topology;
256:   PetscFunctionReturn(PETSC_SUCCESS);
257: }

259: /*@
260:   DMForestSetBaseDM - During the pre-setup phase, set the DM that defines the base mesh of a DMForest forest.  The
261:   forest will be hierarchically refined from the base, and all refinements/coarsenings of the forest will share its
262:   base.  In general, two forest must share a base to be comparable, to do things like construct interpolators.

264:   Logically collective on dm

266:   Input Parameters:
267: + dm - the forest
268: - base - the base DM of the forest

270:   Notes:
271:     Currently the base DM must be a DMPLEX

273:   Level: intermediate

275: .seealso: `DMForestGetBaseDM()`
276: @*/
277: PetscErrorCode DMForestSetBaseDM(DM dm, DM base)
278: {
279:   DM_Forest *forest = (DM_Forest *)dm->data;
280:   PetscInt   dim, dimEmbed;

282:   PetscFunctionBegin;
284:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the base after setup");
285:   PetscCall(PetscObjectReference((PetscObject)base));
286:   PetscCall(DMDestroy(&forest->base));
287:   forest->base = base;
288:   if (base) {
289:     const PetscReal *maxCell, *Lstart, *L;

292:     PetscCall(DMGetDimension(base, &dim));
293:     PetscCall(DMSetDimension(dm, dim));
294:     PetscCall(DMGetCoordinateDim(base, &dimEmbed));
295:     PetscCall(DMSetCoordinateDim(dm, dimEmbed));
296:     PetscCall(DMGetPeriodicity(base, &maxCell, &Lstart, &L));
297:     PetscCall(DMSetPeriodicity(dm, maxCell, Lstart, L));
298:   } else PetscCall(DMSetPeriodicity(dm, NULL, NULL, NULL));
299:   PetscFunctionReturn(PETSC_SUCCESS);
300: }

302: /*@
303:   DMForestGetBaseDM - Get the base DM of a DMForest forest.  The forest will be hierarchically refined from the base,
304:   and all refinements/coarsenings of the forest will share its base.  In general, two forest must share a base to be
305:   comparable, to do things like construct interpolators.

307:   Not collective

309:   Input Parameter:
310: . dm - the forest

312:   Output Parameter:
313: . base - the base DM of the forest

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

318:   Level: intermediate

320: .seealso: `DMForestSetBaseDM()`
321: @*/
322: PetscErrorCode DMForestGetBaseDM(DM dm, DM *base)
323: {
324:   DM_Forest *forest = (DM_Forest *)dm->data;

326:   PetscFunctionBegin;
329:   *base = forest->base;
330:   PetscFunctionReturn(PETSC_SUCCESS);
331: }

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

337:   PetscFunctionBegin;
339:   forest->mapcoordinates    = func;
340:   forest->mapcoordinatesctx = ctx;
341:   PetscFunctionReturn(PETSC_SUCCESS);
342: }

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

348:   PetscFunctionBegin;
350:   if (func) *func = forest->mapcoordinates;
351:   if (ctx) *((void **)ctx) = forest->mapcoordinatesctx;
352:   PetscFunctionReturn(PETSC_SUCCESS);
353: }

355: /*@
356:   DMForestSetAdaptivityForest - During the pre-setup phase, set the forest from which the current forest will be
357:   adapted (e.g., the current forest will be refined/coarsened/repartitioned from it) im DMSetUp().  Usually not needed
358:   by users directly: DMForestTemplate() constructs a new forest to be adapted from an old forest and calls this
359:   routine.

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

364:   Logically collective on dm

366:   Input Parameters:
367: + dm - the new forest, which will be constructed from adapt
368: - adapt - the old forest

370:   Level: intermediate

372: .seealso: `DMForestGetAdaptivityForest()`, `DMForestSetAdaptivityPurpose()`
373: @*/
374: PetscErrorCode DMForestSetAdaptivityForest(DM dm, DM adapt)
375: {
376:   DM_Forest *forest, *adaptForest, *oldAdaptForest;
377:   DM         oldAdapt;
378:   PetscBool  isForest;

380:   PetscFunctionBegin;
383:   PetscCall(DMIsForest(dm, &isForest));
384:   if (!isForest) PetscFunctionReturn(PETSC_SUCCESS);
385:   PetscCheck(adapt == NULL || !dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the adaptation forest after setup");
386:   forest = (DM_Forest *)dm->data;
387:   PetscCall(DMForestGetAdaptivityForest(dm, &oldAdapt));
388:   adaptForest    = (DM_Forest *)(adapt ? adapt->data : NULL);
389:   oldAdaptForest = (DM_Forest *)(oldAdapt ? oldAdapt->data : NULL);
390:   if (adaptForest != oldAdaptForest) {
391:     PetscCall(PetscSFDestroy(&forest->preCoarseToFine));
392:     PetscCall(PetscSFDestroy(&forest->coarseToPreFine));
393:     if (forest->clearadaptivityforest) PetscCall((*forest->clearadaptivityforest)(dm));
394:   }
395:   switch (forest->adaptPurpose) {
396:   case DM_ADAPT_DETERMINE:
397:     PetscCall(PetscObjectReference((PetscObject)adapt));
398:     PetscCall(DMDestroy(&(forest->adapt)));
399:     forest->adapt = adapt;
400:     break;
401:   case DM_ADAPT_REFINE:
402:     PetscCall(DMSetCoarseDM(dm, adapt));
403:     break;
404:   case DM_ADAPT_COARSEN:
405:   case DM_ADAPT_COARSEN_LAST:
406:     PetscCall(DMSetFineDM(dm, adapt));
407:     break;
408:   default:
409:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "invalid adaptivity purpose");
410:   }
411:   PetscFunctionReturn(PETSC_SUCCESS);
412: }

414: /*@
415:   DMForestGetAdaptivityForest - Get the forest from which the current forest is adapted.

417:   Not collective

419:   Input Parameter:
420: . dm - the forest

422:   Output Parameter:
423: . adapt - the forest from which dm is/was adapted

425:   Level: intermediate

427: .seealso: `DMForestSetAdaptivityForest()`, `DMForestSetAdaptivityPurpose()`
428: @*/
429: PetscErrorCode DMForestGetAdaptivityForest(DM dm, DM *adapt)
430: {
431:   DM_Forest *forest;

433:   PetscFunctionBegin;
435:   forest = (DM_Forest *)dm->data;
436:   switch (forest->adaptPurpose) {
437:   case DM_ADAPT_DETERMINE:
438:     *adapt = forest->adapt;
439:     break;
440:   case DM_ADAPT_REFINE:
441:     PetscCall(DMGetCoarseDM(dm, adapt));
442:     break;
443:   case DM_ADAPT_COARSEN:
444:   case DM_ADAPT_COARSEN_LAST:
445:     PetscCall(DMGetFineDM(dm, adapt));
446:     break;
447:   default:
448:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "invalid adaptivity purpose");
449:   }
450:   PetscFunctionReturn(PETSC_SUCCESS);
451: }

453: /*@
454:   DMForestSetAdaptivityPurpose - During the pre-setup phase, set whether the current DM is being adapted from its
455:   source (set with DMForestSetAdaptivityForest()) for the purpose of refinement (DM_ADAPT_REFINE), coarsening
456:   (DM_ADAPT_COARSEN), or undefined (DM_ADAPT_DETERMINE).  This only matters for the purposes of reference counting:
457:   during DMDestroy(), cyclic references can be found between DMs only if the cyclic reference is due to a fine/coarse
458:   relationship (see DMSetFineDM()/DMSetCoarseDM()).  If the purpose is not refinement or coarsening, and the user does
459:   not maintain a reference to the post-adaptation forest (i.e., the one created by DMForestTemplate()), then this can
460:   cause a memory leak.  This method is used by subtypes of DMForest when automatically constructing mesh hierarchies.

462:   Logically collective on dm

464:   Input Parameters:
465: + dm - the forest
466: - purpose - the adaptivity purpose

468:   Level: advanced

470: .seealso: `DMForestTemplate()`, `DMForestSetAdaptivityForest()`, `DMForestGetAdaptivityForest()`, `DMAdaptFlag`
471: @*/
472: PetscErrorCode DMForestSetAdaptivityPurpose(DM dm, DMAdaptFlag purpose)
473: {
474:   DM_Forest *forest;

476:   PetscFunctionBegin;
477:   forest = (DM_Forest *)dm->data;
478:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the adaptation forest after setup");
479:   if (purpose != forest->adaptPurpose) {
480:     DM adapt;

482:     PetscCall(DMForestGetAdaptivityForest(dm, &adapt));
483:     PetscCall(PetscObjectReference((PetscObject)adapt));
484:     PetscCall(DMForestSetAdaptivityForest(dm, NULL));

486:     forest->adaptPurpose = purpose;

488:     PetscCall(DMForestSetAdaptivityForest(dm, adapt));
489:     PetscCall(DMDestroy(&adapt));
490:   }
491:   PetscFunctionReturn(PETSC_SUCCESS);
492: }

494: /*@
495:   DMForestGetAdaptivityPurpose - Get whether the current DM is being adapted from its source (set with
496:   DMForestSetAdaptivityForest()) for the purpose of refinement (DM_ADAPT_REFINE), coarsening (DM_ADAPT_COARSEN),
497:   coarsening only the last level (DM_ADAPT_COARSEN_LAST) or undefined (DM_ADAPT_DETERMINE).
498:   This only matters for the purposes of reference counting: during DMDestroy(), cyclic
499:   references can be found between DMs only if the cyclic reference is due to a fine/coarse relationship (see
500:   DMSetFineDM()/DMSetCoarseDM()).  If the purpose is not refinement or coarsening, and the user does not maintain a
501:   reference to the post-adaptation forest (i.e., the one created by DMForestTemplate()), then this can cause a memory
502:   leak.  This method is used by subtypes of DMForest when automatically constructing mesh hierarchies.

504:   Not collective

506:   Input Parameter:
507: . dm - the forest

509:   Output Parameter:
510: . purpose - the adaptivity purpose

512:   Level: advanced

514: .seealso: `DMForestTemplate()`, `DMForestSetAdaptivityForest()`, `DMForestGetAdaptivityForest()`, `DMAdaptFlag`
515: @*/
516: PetscErrorCode DMForestGetAdaptivityPurpose(DM dm, DMAdaptFlag *purpose)
517: {
518:   DM_Forest *forest;

520:   PetscFunctionBegin;
521:   forest   = (DM_Forest *)dm->data;
522:   *purpose = forest->adaptPurpose;
523:   PetscFunctionReturn(PETSC_SUCCESS);
524: }

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

530:   Logically collective on dm

532:   Input Parameters:
533: + dm - the forest
534: - adjDim - default 0 (i.e., vertices determine adjacency)

536:   Level: intermediate

538: .seealso: `DMForestGetAdjacencyDimension()`, `DMForestSetAdjacencyCodimension()`, `DMForestSetPartitionOverlap()`
539: @*/
540: PetscErrorCode DMForestSetAdjacencyDimension(DM dm, PetscInt adjDim)
541: {
542:   PetscInt   dim;
543:   DM_Forest *forest = (DM_Forest *)dm->data;

545:   PetscFunctionBegin;
547:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the adjacency dimension after setup");
548:   PetscCheck(adjDim >= 0, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "adjacency dim cannot be < 0: %" PetscInt_FMT, adjDim);
549:   PetscCall(DMGetDimension(dm, &dim));
550:   PetscCheck(adjDim <= dim, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "adjacency dim cannot be > %" PetscInt_FMT ": %" PetscInt_FMT, dim, adjDim);
551:   forest->adjDim = adjDim;
552:   PetscFunctionReturn(PETSC_SUCCESS);
553: }

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

559:   Logically collective on dm

561:   Input Parameters:
562: + dm - the forest
563: - adjCodim - default isthe dimension of the forest (see DMGetDimension()), since this is the codimension of vertices

565:   Level: intermediate

567: .seealso: `DMForestGetAdjacencyCodimension()`, `DMForestSetAdjacencyDimension()`
568: @*/
569: PetscErrorCode DMForestSetAdjacencyCodimension(DM dm, PetscInt adjCodim)
570: {
571:   PetscInt dim;

573:   PetscFunctionBegin;
575:   PetscCall(DMGetDimension(dm, &dim));
576:   PetscCall(DMForestSetAdjacencyDimension(dm, dim - adjCodim));
577:   PetscFunctionReturn(PETSC_SUCCESS);
578: }

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

584:   Not collective

586:   Input Parameter:
587: . dm - the forest

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

592:   Level: intermediate

594: .seealso: `DMForestSetAdjacencyDimension()`, `DMForestGetAdjacencyCodimension()`, `DMForestSetPartitionOverlap()`
595: @*/
596: PetscErrorCode DMForestGetAdjacencyDimension(DM dm, PetscInt *adjDim)
597: {
598:   DM_Forest *forest = (DM_Forest *)dm->data;

600:   PetscFunctionBegin;
603:   *adjDim = forest->adjDim;
604:   PetscFunctionReturn(PETSC_SUCCESS);
605: }

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

611:   Not collective

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

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

619:   Level: intermediate

621: .seealso: `DMForestSetAdjacencyCodimension()`, `DMForestGetAdjacencyDimension()`
622: @*/
623: PetscErrorCode DMForestGetAdjacencyCodimension(DM dm, PetscInt *adjCodim)
624: {
625:   DM_Forest *forest = (DM_Forest *)dm->data;
626:   PetscInt   dim;

628:   PetscFunctionBegin;
631:   PetscCall(DMGetDimension(dm, &dim));
632:   *adjCodim = dim - forest->adjDim;
633:   PetscFunctionReturn(PETSC_SUCCESS);
634: }

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

641:   Logically collective on dm

643:   Input Parameters:
644: + dm - the forest
645: - overlap - default 0

647:   Level: intermediate

649: .seealso: `DMForestGetPartitionOverlap()`, `DMForestSetAdjacencyDimension()`, `DMForestSetAdjacencyCodimension()`
650: @*/
651: PetscErrorCode DMForestSetPartitionOverlap(DM dm, PetscInt overlap)
652: {
653:   DM_Forest *forest = (DM_Forest *)dm->data;

655:   PetscFunctionBegin;
657:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the overlap after setup");
658:   PetscCheck(overlap >= 0, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "overlap cannot be < 0: %" PetscInt_FMT, overlap);
659:   forest->overlap = overlap;
660:   PetscFunctionReturn(PETSC_SUCCESS);
661: }

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

667:   Not collective

669:   Input Parameter:
670: . dm - the forest

672:   Output Parameter:
673: . overlap - default 0

675:   Level: intermediate

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

683:   PetscFunctionBegin;
686:   *overlap = forest->overlap;
687:   PetscFunctionReturn(PETSC_SUCCESS);
688: }

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

695:   Logically collective on dm

697:   Input Parameters:
698: + dm - the forest
699: - minRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)

701:   Level: intermediate

703: .seealso: `DMForestGetMinimumRefinement()`, `DMForestSetMaximumRefinement()`, `DMForestSetInitialRefinement()`, `DMForestGetBaseDM()`, `DMForestGetAdaptivityForest()`
704: @*/
705: PetscErrorCode DMForestSetMinimumRefinement(DM dm, PetscInt minRefinement)
706: {
707:   DM_Forest *forest = (DM_Forest *)dm->data;

709:   PetscFunctionBegin;
711:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the minimum refinement after setup");
712:   forest->minRefinement = minRefinement;
713:   PetscFunctionReturn(PETSC_SUCCESS);
714: }

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

721:   Not collective

723:   Input Parameter:
724: . dm - the forest

726:   Output Parameter:
727: . minRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)

729:   Level: intermediate

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

737:   PetscFunctionBegin;
740:   *minRefinement = forest->minRefinement;
741:   PetscFunctionReturn(PETSC_SUCCESS);
742: }

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

748:   Logically collective on dm

750:   Input Parameters:
751: + dm - the forest
752: - initefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)

754:   Level: intermediate

756: .seealso: `DMForestSetMinimumRefinement()`, `DMForestSetMaximumRefinement()`, `DMForestGetBaseDM()`
757: @*/
758: PetscErrorCode DMForestSetInitialRefinement(DM dm, PetscInt initRefinement)
759: {
760:   DM_Forest *forest = (DM_Forest *)dm->data;

762:   PetscFunctionBegin;
764:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the initial refinement after setup");
765:   forest->initRefinement = initRefinement;
766:   PetscFunctionReturn(PETSC_SUCCESS);
767: }

769: /*@
770:   DMForestGetInitialRefinement - Get the initial level of refinement (relative to the base DM, see
771:   DMForestGetBaseDM()) allowed in the forest.

773:   Not collective

775:   Input Parameter:
776: . dm - the forest

778:   Output Parameter:
779: . initRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)

781:   Level: intermediate

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

789:   PetscFunctionBegin;
792:   *initRefinement = forest->initRefinement;
793:   PetscFunctionReturn(PETSC_SUCCESS);
794: }

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

801:   Logically collective on dm

803:   Input Parameters:
804: + dm - the forest
805: - maxRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)

807:   Level: intermediate

809: .seealso: `DMForestGetMinimumRefinement()`, `DMForestSetMaximumRefinement()`, `DMForestSetInitialRefinement()`, `DMForestGetBaseDM()`, `DMForestGetAdaptivityDM()`
810: @*/
811: PetscErrorCode DMForestSetMaximumRefinement(DM dm, PetscInt maxRefinement)
812: {
813:   DM_Forest *forest = (DM_Forest *)dm->data;

815:   PetscFunctionBegin;
817:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the maximum refinement after setup");
818:   forest->maxRefinement = maxRefinement;
819:   PetscFunctionReturn(PETSC_SUCCESS);
820: }

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

827:   Not collective

829:   Input Parameter:
830: . dm - the forest

832:   Output Parameter:
833: . maxRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)

835:   Level: intermediate

837: .seealso: `DMForestSetMaximumRefinement()`, `DMForestGetMinimumRefinement()`, `DMForestGetInitialRefinement()`, `DMForestGetBaseDM()`, `DMForestGetAdaptivityForest()`
838: @*/
839: PetscErrorCode DMForestGetMaximumRefinement(DM dm, PetscInt *maxRefinement)
840: {
841:   DM_Forest *forest = (DM_Forest *)dm->data;

843:   PetscFunctionBegin;
846:   *maxRefinement = forest->maxRefinement;
847:   PetscFunctionReturn(PETSC_SUCCESS);
848: }

850: /*@C
851:   DMForestSetAdaptivityStrategy - During the pre-setup phase, set the strategy for combining adaptivity labels from multiple processes.
852:   Subtypes of DMForest may define their own strategies.  Two default strategies are DMFORESTADAPTALL, which indicates that all processes must agree
853:   for a refinement/coarsening flag to be valid, and DMFORESTADAPTANY, which indicates that only one process needs to
854:   specify refinement/coarsening.

856:   Logically collective on dm

858:   Input Parameters:
859: + dm - the forest
860: - adaptStrategy - default DMFORESTADAPTALL

862:   Level: advanced

864: .seealso: `DMForestGetAdaptivityStrategy()`
865: @*/
866: PetscErrorCode DMForestSetAdaptivityStrategy(DM dm, DMForestAdaptivityStrategy adaptStrategy)
867: {
868:   DM_Forest *forest = (DM_Forest *)dm->data;

870:   PetscFunctionBegin;
872:   PetscCall(PetscFree(forest->adaptStrategy));
873:   PetscCall(PetscStrallocpy((const char *)adaptStrategy, (char **)&forest->adaptStrategy));
874:   PetscFunctionReturn(PETSC_SUCCESS);
875: }

877: /*@C
878:   DMForestSetAdaptivityStrategy - Get the strategy for combining adaptivity labels from multiple processes.  Subtypes
879:   of DMForest may define their own strategies.  Two default strategies are DMFORESTADAPTALL, which indicates that all
880:   processes must agree for a refinement/coarsening flag to be valid, and DMFORESTADAPTANY, which indicates that only
881:   one process needs to specify refinement/coarsening.

883:   Not collective

885:   Input Parameter:
886: . dm - the forest

888:   Output Parameter:
889: . adaptStrategy - the adaptivity strategy (default DMFORESTADAPTALL)

891:   Level: advanced

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

899:   PetscFunctionBegin;
902:   *adaptStrategy = forest->adaptStrategy;
903:   PetscFunctionReturn(PETSC_SUCCESS);
904: }

906: /*@
907:   DMForestGetAdaptivitySuccess - Return whether the requested adaptation (refinement, coarsening, repartitioning,
908:   etc.) was successful.  PETSC_FALSE indicates that the post-adaptation forest is the same as the pre-adpatation
909:   forest.  A requested adaptation may have been unsuccessful if, for example, the requested refinement would have
910:   exceeded the maximum refinement level.

912:   Collective on dm

914:   Input Parameter:

916: . dm - the post-adaptation forest

918:   Output Parameter:

920: . success - PETSC_TRUE if the post-adaptation forest is different from the pre-adaptation forest.

922:   Level: intermediate

924: .see
925: @*/
926: PetscErrorCode DMForestGetAdaptivitySuccess(DM dm, PetscBool *success)
927: {
928:   DM_Forest *forest;

930:   PetscFunctionBegin;
932:   PetscCheck(dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "DMSetUp() has not been called yet.");
933:   forest = (DM_Forest *)dm->data;
934:   PetscCall((forest->getadaptivitysuccess)(dm, success));
935:   PetscFunctionReturn(PETSC_SUCCESS);
936: }

938: /*@
939:   DMForestSetComputeAdaptivitySF - During the pre-setup phase, set whether transfer PetscSFs should be computed
940:   relating the cells of the pre-adaptation forest to the post-adaptiation forest.  After DMSetUp() is called, these transfer PetscSFs can be accessed with DMForestGetAdaptivitySF().

942:   Logically collective on dm

944:   Input Parameters:
945: + dm - the post-adaptation forest
946: - computeSF - default PETSC_TRUE

948:   Level: advanced

950: .seealso: `DMForestGetComputeAdaptivitySF()`, `DMForestGetAdaptivitySF()`
951: @*/
952: PetscErrorCode DMForestSetComputeAdaptivitySF(DM dm, PetscBool computeSF)
953: {
954:   DM_Forest *forest;

956:   PetscFunctionBegin;
958:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot compute adaptivity PetscSFs after setup is called");
959:   forest                 = (DM_Forest *)dm->data;
960:   forest->computeAdaptSF = computeSF;
961:   PetscFunctionReturn(PETSC_SUCCESS);
962: }

964: PetscErrorCode DMForestTransferVec(DM dmIn, Vec vecIn, DM dmOut, Vec vecOut, PetscBool useBCs, PetscReal time)
965: {
966:   DM_Forest *forest;

968:   PetscFunctionBegin;
973:   forest = (DM_Forest *)dmIn->data;
974:   PetscCheck(forest->transfervec, PetscObjectComm((PetscObject)dmIn), PETSC_ERR_SUP, "DMForestTransferVec() not implemented");
975:   PetscCall((forest->transfervec)(dmIn, vecIn, dmOut, vecOut, useBCs, time));
976:   PetscFunctionReturn(PETSC_SUCCESS);
977: }

979: PetscErrorCode DMForestTransferVecFromBase(DM dm, Vec vecIn, Vec vecOut)
980: {
981:   DM_Forest *forest;

983:   PetscFunctionBegin;
987:   forest = (DM_Forest *)dm->data;
988:   PetscCheck(forest->transfervecfrombase, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "DMForestTransferVecFromBase() not implemented");
989:   PetscCall((forest->transfervecfrombase)(dm, vecIn, vecOut));
990:   PetscFunctionReturn(PETSC_SUCCESS);
991: }

993: /*@
994:   DMForestGetComputeAdaptivitySF - Get whether transfer PetscSFs should be computed relating the cells of the
995:   pre-adaptation forest to the post-adaptiation forest.  After DMSetUp() is called, these transfer PetscSFs can be
996:   accessed with DMForestGetAdaptivitySF().

998:   Not collective

1000:   Input Parameter:
1001: . dm - the post-adaptation forest

1003:   Output Parameter:
1004: . computeSF - default PETSC_TRUE

1006:   Level: advanced

1008: .seealso: `DMForestSetComputeAdaptivitySF()`, `DMForestGetAdaptivitySF()`
1009: @*/
1010: PetscErrorCode DMForestGetComputeAdaptivitySF(DM dm, PetscBool *computeSF)
1011: {
1012:   DM_Forest *forest;

1014:   PetscFunctionBegin;
1016:   forest     = (DM_Forest *)dm->data;
1017:   *computeSF = forest->computeAdaptSF;
1018:   PetscFunctionReturn(PETSC_SUCCESS);
1019: }

1021: /*@
1022:   DMForestGetAdaptivitySF - Get PetscSFs that relate the pre-adaptation forest to the post-adaptation forest.
1023:   Adaptation can be any combination of refinement, coarsening, repartition, and change of overlap, so there may be
1024:   some cells of the pre-adaptation that are parents of post-adaptation cells, and vice versa.  Therefore there are two
1025:   PetscSFs: one that relates pre-adaptation coarse cells to post-adaptation fine cells, and one that relates
1026:   pre-adaptation fine cells to post-adaptation coarse cells.

1028:   Not collective

1030:   Input Parameter:
1031:   dm - the post-adaptation forest

1033:   Output Parameter:
1034:   preCoarseToFine - pre-adaptation coarse cells to post-adaptation fine cells: BCast goes from pre- to post-
1035:   coarseToPreFine - post-adaptation coarse cells to pre-adaptation fine cells: BCast goes from post- to pre-

1037:   Level: advanced

1039: .seealso: `DMForestGetComputeAdaptivitySF()`, `DMForestSetComputeAdaptivitySF()`
1040: @*/
1041: PetscErrorCode DMForestGetAdaptivitySF(DM dm, PetscSF *preCoarseToFine, PetscSF *coarseToPreFine)
1042: {
1043:   DM_Forest *forest;

1045:   PetscFunctionBegin;
1047:   PetscCall(DMSetUp(dm));
1048:   forest = (DM_Forest *)dm->data;
1049:   if (preCoarseToFine) *preCoarseToFine = forest->preCoarseToFine;
1050:   if (coarseToPreFine) *coarseToPreFine = forest->coarseToPreFine;
1051:   PetscFunctionReturn(PETSC_SUCCESS);
1052: }

1054: /*@
1055:   DMForestSetGradeFactor - During the pre-setup phase, set the desired amount of grading in the mesh, e.g. give 2 to
1056:   indicate that the diameter of neighboring cells should differ by at most a factor of 2.  Subtypes of DMForest may
1057:   only support one particular choice of grading factor.

1059:   Logically collective on dm

1061:   Input Parameters:
1062: + dm - the forest
1063: - grade - the grading factor

1065:   Level: advanced

1067: .seealso: `DMForestGetGradeFactor()`
1068: @*/
1069: PetscErrorCode DMForestSetGradeFactor(DM dm, PetscInt grade)
1070: {
1071:   DM_Forest *forest = (DM_Forest *)dm->data;

1073:   PetscFunctionBegin;
1075:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the grade factor after setup");
1076:   forest->gradeFactor = grade;
1077:   PetscFunctionReturn(PETSC_SUCCESS);
1078: }

1080: /*@
1081:   DMForestGetGradeFactor - Get the desired amount of grading in the mesh, e.g. give 2 to indicate that the diameter of
1082:   neighboring cells should differ by at most a factor of 2.  Subtypes of DMForest may only support one particular
1083:   choice of grading factor.

1085:   Not collective

1087:   Input Parameter:
1088: . dm - the forest

1090:   Output Parameter:
1091: . grade - the grading factor

1093:   Level: advanced

1095: .seealso: `DMForestSetGradeFactor()`
1096: @*/
1097: PetscErrorCode DMForestGetGradeFactor(DM dm, PetscInt *grade)
1098: {
1099:   DM_Forest *forest = (DM_Forest *)dm->data;

1101:   PetscFunctionBegin;
1104:   *grade = forest->gradeFactor;
1105:   PetscFunctionReturn(PETSC_SUCCESS);
1106: }

1108: /*@
1109:   DMForestSetCellWeightFactor - During the pre-setup phase, set the factor by which the level of refinement changes
1110:   the cell weight (see DMForestSetCellWeights()) when calculating partitions.  The final weight of a cell will be
1111:   (cellWeight) * (weightFactor^refinementLevel).  A factor of 1 indicates that the weight of a cell does not depend on
1112:   its level; a factor of 2, for example, might be appropriate for sub-cycling time-stepping methods, when the
1113:   computation associated with a cell is multiplied by a factor of 2 for each additional level of refinement.

1115:   Logically collective on dm

1117:   Input Parameters:
1118: + dm - the forest
1119: - weightsFactors - default 1.

1121:   Level: advanced

1123: .seealso: `DMForestGetCellWeightFactor()`, `DMForestSetCellWeights()`
1124: @*/
1125: PetscErrorCode DMForestSetCellWeightFactor(DM dm, PetscReal weightsFactor)
1126: {
1127:   DM_Forest *forest = (DM_Forest *)dm->data;

1129:   PetscFunctionBegin;
1131:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the weights factor after setup");
1132:   forest->weightsFactor = weightsFactor;
1133:   PetscFunctionReturn(PETSC_SUCCESS);
1134: }

1136: /*@
1137:   DMForestGetCellWeightFactor - Get the factor by which the level of refinement changes the cell weight (see
1138:   DMForestSetCellWeights()) when calculating partitions.  The final weight of a cell will be (cellWeight) *
1139:   (weightFactor^refinementLevel).  A factor of 1 indicates that the weight of a cell does not depend on its level; a
1140:   factor of 2, for example, might be appropriate for sub-cycling time-stepping methods, when the computation
1141:   associated with a cell is multiplied by a factor of 2 for each additional level of refinement.

1143:   Not collective

1145:   Input Parameter:
1146: . dm - the forest

1148:   Output Parameter:
1149: . weightsFactors - default 1.

1151:   Level: advanced

1153: .seealso: `DMForestSetCellWeightFactor()`, `DMForestSetCellWeights()`
1154: @*/
1155: PetscErrorCode DMForestGetCellWeightFactor(DM dm, PetscReal *weightsFactor)
1156: {
1157:   DM_Forest *forest = (DM_Forest *)dm->data;

1159:   PetscFunctionBegin;
1162:   *weightsFactor = forest->weightsFactor;
1163:   PetscFunctionReturn(PETSC_SUCCESS);
1164: }

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

1169:   Not collective

1171:   Input Parameter:
1172: . dm - the forest

1174:   Output Parameters:
1175: + cStart - the first cell on this process
1176: - cEnd - one after the final cell on this process

1178:   Level: intermediate

1180: .seealso: `DMForestGetCellSF()`
1181: @*/
1182: PetscErrorCode DMForestGetCellChart(DM dm, PetscInt *cStart, PetscInt *cEnd)
1183: {
1184:   DM_Forest *forest = (DM_Forest *)dm->data;

1186:   PetscFunctionBegin;
1190:   if (((forest->cStart == PETSC_DETERMINE) || (forest->cEnd == PETSC_DETERMINE)) && forest->createcellchart) PetscCall(forest->createcellchart(dm, &forest->cStart, &forest->cEnd));
1191:   *cStart = forest->cStart;
1192:   *cEnd   = forest->cEnd;
1193:   PetscFunctionReturn(PETSC_SUCCESS);
1194: }

1196: /*@
1197:   DMForestGetCellSF - After the setup phase, get the PetscSF for overlapping cells between processes

1199:   Not collective

1201:   Input Parameter:
1202: . dm - the forest

1204:   Output Parameter:
1205: . cellSF - the PetscSF

1207:   Level: intermediate

1209: .seealso: `DMForestGetCellChart()`
1210: @*/
1211: PetscErrorCode DMForestGetCellSF(DM dm, PetscSF *cellSF)
1212: {
1213:   DM_Forest *forest = (DM_Forest *)dm->data;

1215:   PetscFunctionBegin;
1218:   if ((!forest->cellSF) && forest->createcellsf) PetscCall(forest->createcellsf(dm, &forest->cellSF));
1219:   *cellSF = forest->cellSF;
1220:   PetscFunctionReturn(PETSC_SUCCESS);
1221: }

1223: /*@C
1224:   DMForestSetAdaptivityLabel - During the pre-setup phase, set the label of the pre-adaptation forest (see
1225:   DMForestGetAdaptivityForest()) that holds the adaptation flags (refinement, coarsening, or some combination).  The
1226:   interpretation of the label values is up to the subtype of DMForest, but DM_ADAPT_DETERMINE, DM_ADAPT_KEEP,
1227:   DM_ADAPT_REFINE, and DM_ADAPT_COARSEN have been reserved as choices that should be accepted by all subtypes.

1229:   Logically collective on dm

1231:   Input Parameters:
1232: - dm - the forest
1233: + adaptLabel - the label in the pre-adaptation forest

1235:   Level: intermediate

1237: .seealso `DMForestGetAdaptivityLabel()`
1238: @*/
1239: PetscErrorCode DMForestSetAdaptivityLabel(DM dm, DMLabel adaptLabel)
1240: {
1241:   DM_Forest *forest = (DM_Forest *)dm->data;

1243:   PetscFunctionBegin;
1246:   PetscCall(PetscObjectReference((PetscObject)adaptLabel));
1247:   PetscCall(DMLabelDestroy(&forest->adaptLabel));
1248:   forest->adaptLabel = adaptLabel;
1249:   PetscFunctionReturn(PETSC_SUCCESS);
1250: }

1252: /*@C
1253:   DMForestGetAdaptivityLabel - Get the label of the pre-adaptation forest (see DMForestGetAdaptivityForest()) that
1254:   holds the adaptation flags (refinement, coarsening, or some combination).  The interpretation of the label values is
1255:   up to the subtype of DMForest, but DM_ADAPT_DETERMINE, DM_ADAPT_KEEP, DM_ADAPT_REFINE, and DM_ADAPT_COARSEN have
1256:   been reserved as choices that should be accepted by all subtypes.

1258:   Not collective

1260:   Input Parameter:
1261: . dm - the forest

1263:   Output Parameter:
1264: . adaptLabel - the name of the label in the pre-adaptation forest

1266:   Level: intermediate

1268: .seealso `DMForestSetAdaptivityLabel()`
1269: @*/
1270: PetscErrorCode DMForestGetAdaptivityLabel(DM dm, DMLabel *adaptLabel)
1271: {
1272:   DM_Forest *forest = (DM_Forest *)dm->data;

1274:   PetscFunctionBegin;
1276:   *adaptLabel = forest->adaptLabel;
1277:   PetscFunctionReturn(PETSC_SUCCESS);
1278: }

1280: /*@
1281:   DMForestSetCellWeights - Set the weights assigned to each of the cells (see DMForestGetCellChart()) of the current
1282:   process: weights are used to determine parallel partitioning.  Partitions will be created so that each process's
1283:   ratio of weight to capacity (see DMForestSetWeightCapacity()) is roughly equal. If NULL, each cell receives a weight
1284:   of 1.

1286:   Logically collective on dm

1288:   Input Parameters:
1289: + dm - the forest
1290: . weights - the array of weights for all cells, or NULL to indicate each cell has weight 1.
1291: - copyMode - how weights should reference weights

1293:   Level: advanced

1295: .seealso: `DMForestGetCellWeights()`, `DMForestSetWeightCapacity()`
1296: @*/
1297: PetscErrorCode DMForestSetCellWeights(DM dm, PetscReal weights[], PetscCopyMode copyMode)
1298: {
1299:   DM_Forest *forest = (DM_Forest *)dm->data;
1300:   PetscInt   cStart, cEnd;

1302:   PetscFunctionBegin;
1304:   PetscCall(DMForestGetCellChart(dm, &cStart, &cEnd));
1305:   PetscCheck(cEnd >= cStart, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "cell chart [%" PetscInt_FMT ",%" PetscInt_FMT ") is not valid", cStart, cEnd);
1306:   if (copyMode == PETSC_COPY_VALUES) {
1307:     if (forest->cellWeightsCopyMode != PETSC_OWN_POINTER || forest->cellWeights == weights) PetscCall(PetscMalloc1(cEnd - cStart, &forest->cellWeights));
1308:     PetscCall(PetscArraycpy(forest->cellWeights, weights, cEnd - cStart));
1309:     forest->cellWeightsCopyMode = PETSC_OWN_POINTER;
1310:     PetscFunctionReturn(PETSC_SUCCESS);
1311:   }
1312:   if (forest->cellWeightsCopyMode == PETSC_OWN_POINTER) PetscCall(PetscFree(forest->cellWeights));
1313:   forest->cellWeights         = weights;
1314:   forest->cellWeightsCopyMode = copyMode;
1315:   PetscFunctionReturn(PETSC_SUCCESS);
1316: }

1318: /*@
1319:   DMForestGetCellWeights - Get the weights assigned to each of the cells (see DMForestGetCellChart()) of the current
1320:   process: weights are used to determine parallel partitioning.  Partitions will be created so that each process's
1321:   ratio of weight to capacity (see DMForestSetWeightCapacity()) is roughly equal. If NULL, each cell receives a weight
1322:   of 1.

1324:   Not collective

1326:   Input Parameter:
1327: . dm - the forest

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

1332:   Level: advanced

1334: .seealso: `DMForestSetCellWeights()`, `DMForestSetWeightCapacity()`
1335: @*/
1336: PetscErrorCode DMForestGetCellWeights(DM dm, PetscReal **weights)
1337: {
1338:   DM_Forest *forest = (DM_Forest *)dm->data;

1340:   PetscFunctionBegin;
1343:   *weights = forest->cellWeights;
1344:   PetscFunctionReturn(PETSC_SUCCESS);
1345: }

1347: /*@
1348:   DMForestSetWeightCapacity - During the pre-setup phase, set the capacity of the current process when repartitioning
1349:   a pre-adaptation forest (see DMForestGetAdaptivityForest()).  After partitioning, the ratio of the weight of each
1350:   process's cells to the process's capacity will be roughly equal for all processes.  A capacity of 0 indicates that
1351:   the current process should not have any cells after repartitioning.

1353:   Logically Collective on dm

1355:   Input parameters:
1356: + dm - the forest
1357: - capacity - this process's capacity

1359:   Level: advanced

1361: .seealso `DMForestGetWeightCapacity()`, `DMForestSetCellWeights()`, `DMForestSetCellWeightFactor()`
1362: @*/
1363: PetscErrorCode DMForestSetWeightCapacity(DM dm, PetscReal capacity)
1364: {
1365:   DM_Forest *forest = (DM_Forest *)dm->data;

1367:   PetscFunctionBegin;
1369:   PetscCheck(!dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Cannot change the weight capacity after setup");
1370:   PetscCheck(capacity >= 0., PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "Cannot have negative weight capacity; %g", (double)capacity);
1371:   forest->weightCapacity = capacity;
1372:   PetscFunctionReturn(PETSC_SUCCESS);
1373: }

1375: /*@
1376:   DMForestGetWeightCapacity - Set the capacity of the current process when repartitioning a pre-adaptation forest (see
1377:   DMForestGetAdaptivityForest()).  After partitioning, the ratio of the weight of each process's cells to the
1378:   process's capacity will be roughly equal for all processes.  A capacity of 0 indicates that the current process
1379:   should not have any cells after repartitioning.

1381:   Not collective

1383:   Input parameter:
1384: . dm - the forest

1386:   Output parameter:
1387: . capacity - this process's capacity

1389:   Level: advanced

1391: .seealso `DMForestSetWeightCapacity()`, `DMForestSetCellWeights()`, `DMForestSetCellWeightFactor()`
1392: @*/
1393: PetscErrorCode DMForestGetWeightCapacity(DM dm, PetscReal *capacity)
1394: {
1395:   DM_Forest *forest = (DM_Forest *)dm->data;

1397:   PetscFunctionBegin;
1400:   *capacity = forest->weightCapacity;
1401:   PetscFunctionReturn(PETSC_SUCCESS);
1402: }

1404: PETSC_EXTERN PetscErrorCode DMSetFromOptions_Forest(DM dm, PetscOptionItems *PetscOptionsObject)
1405: {
1406:   PetscBool                  flg, flg1, flg2, flg3, flg4;
1407:   DMForestTopology           oldTopo;
1408:   char                       stringBuffer[256];
1409:   PetscViewer                viewer;
1410:   PetscViewerFormat          format;
1411:   PetscInt                   adjDim, adjCodim, overlap, minRefinement, initRefinement, maxRefinement, grade;
1412:   PetscReal                  weightsFactor;
1413:   DMForestAdaptivityStrategy adaptStrategy;

1415:   PetscFunctionBegin;
1416:   PetscCall(DMForestGetTopology(dm, &oldTopo));
1417:   PetscOptionsHeadBegin(PetscOptionsObject, "DMForest Options");
1418:   PetscCall(PetscOptionsString("-dm_forest_topology", "the topology of the forest's base mesh", "DMForestSetTopology", oldTopo, stringBuffer, sizeof(stringBuffer), &flg1));
1419:   PetscCall(PetscOptionsViewer("-dm_forest_base_dm", "load the base DM from a viewer specification", "DMForestSetBaseDM", &viewer, &format, &flg2));
1420:   PetscCall(PetscOptionsViewer("-dm_forest_coarse_forest", "load the coarse forest from a viewer specification", "DMForestSetCoarseForest", &viewer, &format, &flg3));
1421:   PetscCall(PetscOptionsViewer("-dm_forest_fine_forest", "load the fine forest from a viewer specification", "DMForestSetFineForest", &viewer, &format, &flg4));
1422:   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}");
1423:   if (flg1) {
1424:     PetscCall(DMForestSetTopology(dm, (DMForestTopology)stringBuffer));
1425:     PetscCall(DMForestSetBaseDM(dm, NULL));
1426:     PetscCall(DMForestSetAdaptivityForest(dm, NULL));
1427:   }
1428:   if (flg2) {
1429:     DM base;

1431:     PetscCall(DMCreate(PetscObjectComm((PetscObject)dm), &base));
1432:     PetscCall(PetscViewerPushFormat(viewer, format));
1433:     PetscCall(DMLoad(base, viewer));
1434:     PetscCall(PetscViewerDestroy(&viewer));
1435:     PetscCall(DMForestSetBaseDM(dm, base));
1436:     PetscCall(DMDestroy(&base));
1437:     PetscCall(DMForestSetTopology(dm, NULL));
1438:     PetscCall(DMForestSetAdaptivityForest(dm, NULL));
1439:   }
1440:   if (flg3) {
1441:     DM coarse;

1443:     PetscCall(DMCreate(PetscObjectComm((PetscObject)dm), &coarse));
1444:     PetscCall(PetscViewerPushFormat(viewer, format));
1445:     PetscCall(DMLoad(coarse, viewer));
1446:     PetscCall(PetscViewerDestroy(&viewer));
1447:     PetscCall(DMForestSetAdaptivityForest(dm, coarse));
1448:     PetscCall(DMDestroy(&coarse));
1449:     PetscCall(DMForestSetTopology(dm, NULL));
1450:     PetscCall(DMForestSetBaseDM(dm, NULL));
1451:   }
1452:   if (flg4) {
1453:     DM fine;

1455:     PetscCall(DMCreate(PetscObjectComm((PetscObject)dm), &fine));
1456:     PetscCall(PetscViewerPushFormat(viewer, format));
1457:     PetscCall(DMLoad(fine, viewer));
1458:     PetscCall(PetscViewerDestroy(&viewer));
1459:     PetscCall(DMForestSetAdaptivityForest(dm, fine));
1460:     PetscCall(DMDestroy(&fine));
1461:     PetscCall(DMForestSetTopology(dm, NULL));
1462:     PetscCall(DMForestSetBaseDM(dm, NULL));
1463:   }
1464:   PetscCall(DMForestGetAdjacencyDimension(dm, &adjDim));
1465:   PetscCall(PetscOptionsBoundedInt("-dm_forest_adjacency_dimension", "set the dimension of points that define adjacency in the forest", "DMForestSetAdjacencyDimension", adjDim, &adjDim, &flg, 0));
1466:   if (flg) {
1467:     PetscCall(DMForestSetAdjacencyDimension(dm, adjDim));
1468:   } else {
1469:     PetscCall(DMForestGetAdjacencyCodimension(dm, &adjCodim));
1470:     PetscCall(PetscOptionsBoundedInt("-dm_forest_adjacency_codimension", "set the codimension of points that define adjacency in the forest", "DMForestSetAdjacencyCodimension", adjCodim, &adjCodim, &flg, 1));
1471:     if (flg) PetscCall(DMForestSetAdjacencyCodimension(dm, adjCodim));
1472:   }
1473:   PetscCall(DMForestGetPartitionOverlap(dm, &overlap));
1474:   PetscCall(PetscOptionsBoundedInt("-dm_forest_partition_overlap", "set the degree of partition overlap", "DMForestSetPartitionOverlap", overlap, &overlap, &flg, 0));
1475:   if (flg) PetscCall(DMForestSetPartitionOverlap(dm, overlap));
1476: #if 0
1477:   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));
1478:   if (flg) {
1479:     PetscCall(DMForestSetMinimumRefinement(dm,minRefinement));
1480:     PetscCall(DMForestSetInitialRefinement(dm,minRefinement));
1481:   }
1482:   PetscCall(PetscOptionsBoundedInt("-dm_refine_hierarchy","equivalent to -dm_forest_set_minimum_refinement 0 and -dm_forest_set_initial_refinement",NULL,initRefinement,&initRefinement,&flg,0));
1483:   if (flg) {
1484:     PetscCall(DMForestSetMinimumRefinement(dm,0));
1485:     PetscCall(DMForestSetInitialRefinement(dm,initRefinement));
1486:   }
1487: #endif
1488:   PetscCall(DMForestGetMinimumRefinement(dm, &minRefinement));
1489:   PetscCall(PetscOptionsBoundedInt("-dm_forest_minimum_refinement", "set the minimum level of refinement in the forest", "DMForestSetMinimumRefinement", minRefinement, &minRefinement, &flg, 0));
1490:   if (flg) PetscCall(DMForestSetMinimumRefinement(dm, minRefinement));
1491:   PetscCall(DMForestGetInitialRefinement(dm, &initRefinement));
1492:   PetscCall(PetscOptionsBoundedInt("-dm_forest_initial_refinement", "set the initial level of refinement in the forest", "DMForestSetInitialRefinement", initRefinement, &initRefinement, &flg, 0));
1493:   if (flg) PetscCall(DMForestSetInitialRefinement(dm, initRefinement));
1494:   PetscCall(DMForestGetMaximumRefinement(dm, &maxRefinement));
1495:   PetscCall(PetscOptionsBoundedInt("-dm_forest_maximum_refinement", "set the maximum level of refinement in the forest", "DMForestSetMaximumRefinement", maxRefinement, &maxRefinement, &flg, 0));
1496:   if (flg) PetscCall(DMForestSetMaximumRefinement(dm, maxRefinement));
1497:   PetscCall(DMForestGetAdaptivityStrategy(dm, &adaptStrategy));
1498:   PetscCall(PetscOptionsString("-dm_forest_adaptivity_strategy", "the forest's adaptivity-flag resolution strategy", "DMForestSetAdaptivityStrategy", adaptStrategy, stringBuffer, sizeof(stringBuffer), &flg));
1499:   if (flg) PetscCall(DMForestSetAdaptivityStrategy(dm, (DMForestAdaptivityStrategy)stringBuffer));
1500:   PetscCall(DMForestGetGradeFactor(dm, &grade));
1501:   PetscCall(PetscOptionsBoundedInt("-dm_forest_grade_factor", "grade factor between neighboring cells", "DMForestSetGradeFactor", grade, &grade, &flg, 0));
1502:   if (flg) PetscCall(DMForestSetGradeFactor(dm, grade));
1503:   PetscCall(DMForestGetCellWeightFactor(dm, &weightsFactor));
1504:   PetscCall(PetscOptionsReal("-dm_forest_cell_weight_factor", "multiplying weight factor for cell refinement", "DMForestSetCellWeightFactor", weightsFactor, &weightsFactor, &flg));
1505:   if (flg) PetscCall(DMForestSetCellWeightFactor(dm, weightsFactor));
1506:   PetscOptionsHeadEnd();
1507:   PetscFunctionReturn(PETSC_SUCCESS);
1508: }

1510: PetscErrorCode DMCreateSubDM_Forest(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1511: {
1512:   PetscFunctionBegin;
1513:   if (subdm) PetscCall(DMClone(dm, subdm));
1514:   PetscCall(DMCreateSectionSubDM(dm, numFields, fields, is, subdm));
1515:   PetscFunctionReturn(PETSC_SUCCESS);
1516: }

1518: PetscErrorCode DMRefine_Forest(DM dm, MPI_Comm comm, DM *dmRefined)
1519: {
1520:   DMLabel refine;
1521:   DM      fineDM;

1523:   PetscFunctionBegin;
1524:   PetscCall(DMGetFineDM(dm, &fineDM));
1525:   if (fineDM) {
1526:     PetscCall(PetscObjectReference((PetscObject)fineDM));
1527:     *dmRefined = fineDM;
1528:     PetscFunctionReturn(PETSC_SUCCESS);
1529:   }
1530:   PetscCall(DMForestTemplate(dm, comm, dmRefined));
1531:   PetscCall(DMGetLabel(dm, "refine", &refine));
1532:   if (!refine) {
1533:     PetscCall(DMLabelCreate(PETSC_COMM_SELF, "refine", &refine));
1534:     PetscCall(DMLabelSetDefaultValue(refine, DM_ADAPT_REFINE));
1535:   } else PetscCall(PetscObjectReference((PetscObject)refine));
1536:   PetscCall(DMForestSetAdaptivityLabel(*dmRefined, refine));
1537:   PetscCall(DMLabelDestroy(&refine));
1538:   PetscFunctionReturn(PETSC_SUCCESS);
1539: }

1541: PetscErrorCode DMCoarsen_Forest(DM dm, MPI_Comm comm, DM *dmCoarsened)
1542: {
1543:   DMLabel coarsen;
1544:   DM      coarseDM;

1546:   PetscFunctionBegin;
1547:   {
1548:     PetscMPIInt mpiComparison;
1549:     MPI_Comm    dmcomm = PetscObjectComm((PetscObject)dm);

1551:     PetscCallMPI(MPI_Comm_compare(comm, dmcomm, &mpiComparison));
1552:     PetscCheck(mpiComparison == MPI_IDENT || mpiComparison == MPI_CONGRUENT, dmcomm, PETSC_ERR_SUP, "No support for different communicators yet");
1553:   }
1554:   PetscCall(DMGetCoarseDM(dm, &coarseDM));
1555:   if (coarseDM) {
1556:     PetscCall(PetscObjectReference((PetscObject)coarseDM));
1557:     *dmCoarsened = coarseDM;
1558:     PetscFunctionReturn(PETSC_SUCCESS);
1559:   }
1560:   PetscCall(DMForestTemplate(dm, comm, dmCoarsened));
1561:   PetscCall(DMForestSetAdaptivityPurpose(*dmCoarsened, DM_ADAPT_COARSEN));
1562:   PetscCall(DMGetLabel(dm, "coarsen", &coarsen));
1563:   if (!coarsen) {
1564:     PetscCall(DMLabelCreate(PETSC_COMM_SELF, "coarsen", &coarsen));
1565:     PetscCall(DMLabelSetDefaultValue(coarsen, DM_ADAPT_COARSEN));
1566:   } else PetscCall(PetscObjectReference((PetscObject)coarsen));
1567:   PetscCall(DMForestSetAdaptivityLabel(*dmCoarsened, coarsen));
1568:   PetscCall(DMLabelDestroy(&coarsen));
1569:   PetscFunctionReturn(PETSC_SUCCESS);
1570: }

1572: PetscErrorCode DMAdaptLabel_Forest(DM dm, PETSC_UNUSED Vec metric, DMLabel label, PETSC_UNUSED DMLabel rgLabel, DM *adaptedDM)
1573: {
1574:   PetscBool success;

1576:   PetscFunctionBegin;
1577:   PetscCall(DMForestTemplate(dm, PetscObjectComm((PetscObject)dm), adaptedDM));
1578:   PetscCall(DMForestSetAdaptivityLabel(*adaptedDM, label));
1579:   PetscCall(DMSetUp(*adaptedDM));
1580:   PetscCall(DMForestGetAdaptivitySuccess(*adaptedDM, &success));
1581:   if (!success) {
1582:     PetscCall(DMDestroy(adaptedDM));
1583:     *adaptedDM = NULL;
1584:   }
1585:   PetscFunctionReturn(PETSC_SUCCESS);
1586: }

1588: static PetscErrorCode DMInitialize_Forest(DM dm)
1589: {
1590:   PetscFunctionBegin;
1591:   PetscCall(PetscMemzero(dm->ops, sizeof(*(dm->ops))));

1593:   dm->ops->clone          = DMClone_Forest;
1594:   dm->ops->setfromoptions = DMSetFromOptions_Forest;
1595:   dm->ops->destroy        = DMDestroy_Forest;
1596:   dm->ops->createsubdm    = DMCreateSubDM_Forest;
1597:   dm->ops->refine         = DMRefine_Forest;
1598:   dm->ops->coarsen        = DMCoarsen_Forest;
1599:   PetscFunctionReturn(PETSC_SUCCESS);
1600: }

1602: /*MC

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

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

1615:   Level: advanced

1617: .seealso: `DMType`, `DMCreate()`, `DMSetType()`, `DMForestGetBaseDM()`, `DMForestSetBaseDM()`, `DMForestTemplate()`, `DMForestSetAdaptivityLabel()`
1618: M*/

1620: PETSC_EXTERN PetscErrorCode DMCreate_Forest(DM dm)
1621: {
1622:   DM_Forest *forest;

1624:   PetscFunctionBegin;
1626:   PetscCall(PetscNew(&forest));
1627:   dm->dim                     = 0;
1628:   dm->data                    = forest;
1629:   forest->refct               = 1;
1630:   forest->data                = NULL;
1631:   forest->topology            = NULL;
1632:   forest->adapt               = NULL;
1633:   forest->base                = NULL;
1634:   forest->adaptPurpose        = DM_ADAPT_DETERMINE;
1635:   forest->adjDim              = PETSC_DEFAULT;
1636:   forest->overlap             = PETSC_DEFAULT;
1637:   forest->minRefinement       = PETSC_DEFAULT;
1638:   forest->maxRefinement       = PETSC_DEFAULT;
1639:   forest->initRefinement      = PETSC_DEFAULT;
1640:   forest->cStart              = PETSC_DETERMINE;
1641:   forest->cEnd                = PETSC_DETERMINE;
1642:   forest->cellSF              = NULL;
1643:   forest->adaptLabel          = NULL;
1644:   forest->gradeFactor         = 2;
1645:   forest->cellWeights         = NULL;
1646:   forest->cellWeightsCopyMode = PETSC_USE_POINTER;
1647:   forest->weightsFactor       = 1.;
1648:   forest->weightCapacity      = 1.;
1649:   PetscCall(DMForestSetAdaptivityStrategy(dm, DMFORESTADAPTALL));
1650:   PetscCall(DMInitialize_Forest(dm));
1651:   PetscFunctionReturn(PETSC_SUCCESS);
1652: }