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: {
12: char *name;
13: DMForestTypeLink next;
14: };
16: DMForestTypeLink DMForestTypeList;
18: static PetscErrorCode DMForestPackageFinalize(void)
19: {
20: DMForestTypeLink oldLink, link = DMForestTypeList;
22: while (link) {
23: oldLink = link;
24: PetscFree(oldLink->name);
25: link = oldLink->next;
26: PetscFree(oldLink);
27: }
28: return 0;
29: }
31: static PetscErrorCode DMForestPackageInitialize(void)
32: {
33: if (DMForestPackageInitialized) return 0;
34: DMForestPackageInitialized = PETSC_TRUE;
36: DMForestRegisterType(DMFOREST);
37: PetscRegisterFinalize(DMForestPackageFinalize);
38: return 0;
39: }
41: /*@C
42: DMForestRegisterType - Registers a DMType as a subtype of DMFOREST (so that DMIsForest() will be correct)
44: Not Collective
46: Input parameter:
47: . name - the name of the type
49: Level: advanced
51: .seealso: DMFOREST, DMIsForest()
52: @*/
53: PetscErrorCode DMForestRegisterType(DMType name)
54: {
55: DMForestTypeLink link;
57: DMForestPackageInitialize();
58: PetscNew(&link);
59: PetscStrallocpy(name,&link->name);
60: link->next = DMForestTypeList;
61: DMForestTypeList = link;
62: return 0;
63: }
65: /*@
66: DMIsForest - Check whether a DM uses the DMFOREST interface for hierarchically-refined meshes
68: Not Collective
70: Input parameter:
71: . dm - the DM object
73: Output parameter:
74: . isForest - whether dm is a subtype of DMFOREST
76: Level: intermediate
78: .seealso: DMFOREST, DMForestRegisterType()
79: @*/
80: PetscErrorCode DMIsForest(DM dm, PetscBool *isForest)
81: {
82: DMForestTypeLink link = DMForestTypeList;
84: while (link) {
85: PetscBool sameType;
86: PetscObjectTypeCompare((PetscObject)dm,link->name,&sameType);
87: if (sameType) {
88: *isForest = PETSC_TRUE;
89: return 0;
90: }
91: link = link->next;
92: }
93: *isForest = PETSC_FALSE;
94: return 0;
95: }
97: /*@
98: DMForestTemplate - Create a new DM that will be adapted from a source DM. The new DM reproduces the configuration
99: of the source, but is not yet setup, so that the user can then define only the ways that the new DM should differ
100: (by, e.g., refinement or repartitioning). The source DM is also set as the adaptivity source DM of the new DM (see
101: DMForestSetAdaptivityForest()).
103: Collective on dm
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: .seealso: DMForestSetAdaptivityForest()
115: @*/
116: PetscErrorCode DMForestTemplate(DM dm, MPI_Comm comm, DM *tdm)
117: {
118: DM_Forest *forest = (DM_Forest*) dm->data;
119: DMType type;
120: DM base;
121: DMForestTopology topology;
122: MatType mtype;
123: PetscInt dim, overlap, ref, factor;
124: DMForestAdaptivityStrategy strat;
125: void *ctx;
126: PetscErrorCode (*map)(DM, PetscInt, PetscInt, const PetscReal[], PetscReal[], void*);
127: void *mapCtx;
130: DMCreate(PetscObjectComm((PetscObject)dm),tdm);
131: DMGetType(dm,&type);
132: DMSetType(*tdm,type);
133: DMForestGetBaseDM(dm,&base);
134: DMForestSetBaseDM(*tdm,base);
135: DMForestGetTopology(dm,&topology);
136: DMForestSetTopology(*tdm,topology);
137: DMForestGetAdjacencyDimension(dm,&dim);
138: DMForestSetAdjacencyDimension(*tdm,dim);
139: DMForestGetPartitionOverlap(dm,&overlap);
140: DMForestSetPartitionOverlap(*tdm,overlap);
141: DMForestGetMinimumRefinement(dm,&ref);
142: DMForestSetMinimumRefinement(*tdm,ref);
143: DMForestGetMaximumRefinement(dm,&ref);
144: DMForestSetMaximumRefinement(*tdm,ref);
145: DMForestGetAdaptivityStrategy(dm,&strat);
146: DMForestSetAdaptivityStrategy(*tdm,strat);
147: DMForestGetGradeFactor(dm,&factor);
148: DMForestSetGradeFactor(*tdm,factor);
149: DMForestGetBaseCoordinateMapping(dm,&map,&mapCtx);
150: DMForestSetBaseCoordinateMapping(*tdm,map,mapCtx);
151: if (forest->ftemplate) {
152: (*forest->ftemplate)(dm, *tdm);
153: }
154: DMForestSetAdaptivityForest(*tdm,dm);
155: DMCopyDisc(dm,*tdm);
156: DMGetApplicationContext(dm,&ctx);
157: DMSetApplicationContext(*tdm,&ctx);
158: {
159: PetscBool isper;
160: const PetscReal *maxCell, *L;
161: const DMBoundaryType *bd;
163: DMGetPeriodicity(dm,&isper,&maxCell,&L,&bd);
164: DMSetPeriodicity(*tdm,isper,maxCell,L,bd);
165: }
166: DMGetMatType(dm,&mtype);
167: DMSetMatType(*tdm,mtype);
168: return 0;
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: forest->refct++;
179: (*newdm)->data = forest;
180: PetscObjectGetType((PetscObject) dm, &type);
181: PetscObjectChangeTypeName((PetscObject) *newdm, type);
182: DMInitialize_Forest(*newdm);
183: return 0;
184: }
186: static PetscErrorCode DMDestroy_Forest(DM dm)
187: {
188: DM_Forest *forest = (DM_Forest*) dm->data;
190: if (--forest->refct > 0) return 0;
191: if (forest->destroy) (*forest->destroy)(dm);
192: PetscSFDestroy(&forest->cellSF);
193: PetscSFDestroy(&forest->preCoarseToFine);
194: PetscSFDestroy(&forest->coarseToPreFine);
195: DMLabelDestroy(&forest->adaptLabel);
196: PetscFree(forest->adaptStrategy);
197: DMDestroy(&forest->base);
198: DMDestroy(&forest->adapt);
199: PetscFree(forest->topology);
200: PetscFree(forest);
201: return 0;
202: }
204: /*@C
205: DMForestSetTopology - Set the topology of a DMForest during the pre-setup phase. The topology is a string (e.g.
206: "cube", "shell") and can be interpreted by subtypes of DMFOREST) to construct the base DM of a forest during
207: DMSetUp().
209: Logically collective on dm
211: Input parameters:
212: + dm - the forest
213: - topology - the topology of the forest
215: Level: intermediate
217: .seealso: DMForestGetTopology(), DMForestSetBaseDM()
218: @*/
219: PetscErrorCode DMForestSetTopology(DM dm, DMForestTopology topology)
220: {
221: DM_Forest *forest = (DM_Forest*) dm->data;
225: PetscFree(forest->topology);
226: PetscStrallocpy((const char*)topology,(char**) &forest->topology);
227: return 0;
228: }
230: /*@C
231: DMForestGetTopology - Get a string describing the topology of a DMForest.
233: Not collective
235: Input parameter:
236: . dm - the forest
238: Output parameter:
239: . topology - the topology of the forest (e.g., 'cube', 'shell')
241: Level: intermediate
243: .seealso: DMForestSetTopology()
244: @*/
245: PetscErrorCode DMForestGetTopology(DM dm, DMForestTopology *topology)
246: {
247: DM_Forest *forest = (DM_Forest*) dm->data;
251: *topology = forest->topology;
252: return 0;
253: }
255: /*@
256: DMForestSetBaseDM - During the pre-setup phase, set the DM that defines the base mesh of a DMForest forest. The
257: forest will be hierarchically refined from the base, and all refinements/coarsenings of the forest will share its
258: base. In general, two forest must share a base to be comparable, to do things like construct interpolators.
260: Logically collective on dm
262: Input Parameters:
263: + dm - the forest
264: - base - the base DM of the forest
266: Notes:
267: Currently the base DM must be a DMPLEX
269: Level: intermediate
271: .seealso: DMForestGetBaseDM()
272: @*/
273: PetscErrorCode DMForestSetBaseDM(DM dm, DM base)
274: {
275: DM_Forest *forest = (DM_Forest*) dm->data;
276: PetscInt dim, dimEmbed;
280: PetscObjectReference((PetscObject)base);
281: DMDestroy(&forest->base);
282: forest->base = base;
283: if (base) {
284: PetscBool isper;
285: const PetscReal *maxCell, *L;
286: const DMBoundaryType *bd;
289: DMGetDimension(base,&dim);
290: DMSetDimension(dm,dim);
291: DMGetCoordinateDim(base,&dimEmbed);
292: DMSetCoordinateDim(dm,dimEmbed);
293: DMGetPeriodicity(base,&isper,&maxCell,&L,&bd);
294: DMSetPeriodicity(dm,isper,maxCell,L,bd);
295: } else {
296: DMSetPeriodicity(dm,PETSC_FALSE,NULL,NULL,NULL);
297: }
298: return 0;
299: }
301: /*@
302: DMForestGetBaseDM - Get the base DM of a DMForest forest. The forest will be hierarchically refined from the base,
303: and all refinements/coarsenings of the forest will share its base. In general, two forest must share a base to be
304: comparable, to do things like construct interpolators.
306: Not collective
308: Input Parameter:
309: . dm - the forest
311: Output Parameter:
312: . base - the base DM of the forest
314: Notes:
315: After DMSetUp(), the base DM will be redundantly distributed across MPI processes
317: Level: intermediate
319: .seealso: DMForestSetBaseDM()
320: @*/
321: PetscErrorCode DMForestGetBaseDM(DM dm, DM *base)
322: {
323: DM_Forest *forest = (DM_Forest*) dm->data;
327: *base = forest->base;
328: return 0;
329: }
331: PetscErrorCode DMForestSetBaseCoordinateMapping(DM dm, PetscErrorCode (*func)(DM,PetscInt,PetscInt,const PetscReal [],PetscReal [],void*),void *ctx)
332: {
333: DM_Forest *forest = (DM_Forest*) dm->data;
336: forest->mapcoordinates = func;
337: forest->mapcoordinatesctx = ctx;
338: return 0;
339: }
341: PetscErrorCode DMForestGetBaseCoordinateMapping(DM dm, PetscErrorCode (**func) (DM,PetscInt,PetscInt,const PetscReal [],PetscReal [],void*),void *ctx)
342: {
343: DM_Forest *forest = (DM_Forest*) dm->data;
346: if (func) *func = forest->mapcoordinates;
347: if (ctx) *((void**) ctx) = forest->mapcoordinatesctx;
348: return 0;
349: }
351: /*@
352: DMForestSetAdaptivityForest - During the pre-setup phase, set the forest from which the current forest will be
353: adapted (e.g., the current forest will be refined/coarsened/repartitioned from it) im DMSetUp(). Usually not needed
354: by users directly: DMForestTemplate() constructs a new forest to be adapted from an old forest and calls this
355: routine.
357: Note that this can be called after setup with adapt = NULL, which will clear all internal data related to the
358: adaptivity forest from dm. This way, repeatedly adapting does not leave stale DM objects in memory.
360: Logically collective on dm
362: Input Parameters:
363: + dm - the new forest, which will be constructed from adapt
364: - adapt - the old forest
366: Level: intermediate
368: .seealso: DMForestGetAdaptivityForest(), DMForestSetAdaptivityPurpose()
369: @*/
370: PetscErrorCode DMForestSetAdaptivityForest(DM dm,DM adapt)
371: {
372: DM_Forest *forest, *adaptForest, *oldAdaptForest;
373: DM oldAdapt;
374: PetscBool isForest;
378: DMIsForest(dm, &isForest);
379: if (!isForest) return 0;
381: forest = (DM_Forest*) dm->data;
382: DMForestGetAdaptivityForest(dm,&oldAdapt);
383: adaptForest = (DM_Forest*) (adapt ? adapt->data : NULL);
384: oldAdaptForest = (DM_Forest*) (oldAdapt ? oldAdapt->data : NULL);
385: if (adaptForest != oldAdaptForest) {
386: PetscSFDestroy(&forest->preCoarseToFine);
387: PetscSFDestroy(&forest->coarseToPreFine);
388: if (forest->clearadaptivityforest) (*forest->clearadaptivityforest)(dm);
389: }
390: switch (forest->adaptPurpose) {
391: case DM_ADAPT_DETERMINE:
392: PetscObjectReference((PetscObject)adapt);
393: DMDestroy(&(forest->adapt));
394: forest->adapt = adapt;
395: break;
396: case DM_ADAPT_REFINE:
397: DMSetCoarseDM(dm,adapt);
398: break;
399: case DM_ADAPT_COARSEN:
400: case DM_ADAPT_COARSEN_LAST:
401: DMSetFineDM(dm,adapt);
402: break;
403: default:
404: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"invalid adaptivity purpose");
405: }
406: return 0;
407: }
409: /*@
410: DMForestGetAdaptivityForest - Get the forest from which the current forest is adapted.
412: Not collective
414: Input Parameter:
415: . dm - the forest
417: Output Parameter:
418: . adapt - the forest from which dm is/was adapted
420: Level: intermediate
422: .seealso: DMForestSetAdaptivityForest(), DMForestSetAdaptivityPurpose()
423: @*/
424: PetscErrorCode DMForestGetAdaptivityForest(DM dm, DM *adapt)
425: {
426: DM_Forest *forest;
429: forest = (DM_Forest*) dm->data;
430: switch (forest->adaptPurpose) {
431: case DM_ADAPT_DETERMINE:
432: *adapt = forest->adapt;
433: break;
434: case DM_ADAPT_REFINE:
435: DMGetCoarseDM(dm,adapt);
436: break;
437: case DM_ADAPT_COARSEN:
438: case DM_ADAPT_COARSEN_LAST:
439: DMGetFineDM(dm,adapt);
440: break;
441: default:
442: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"invalid adaptivity purpose");
443: }
444: return 0;
445: }
447: /*@
448: DMForestSetAdaptivityPurpose - During the pre-setup phase, set whether the current DM is being adapted from its
449: source (set with DMForestSetAdaptivityForest()) for the purpose of refinement (DM_ADAPT_REFINE), coarsening
450: (DM_ADAPT_COARSEN), or undefined (DM_ADAPT_DETERMINE). This only matters for the purposes of reference counting:
451: during DMDestroy(), cyclic references can be found between DMs only if the cyclic reference is due to a fine/coarse
452: relationship (see DMSetFineDM()/DMSetCoarseDM()). If the purpose is not refinement or coarsening, and the user does
453: not maintain a reference to the post-adaptation forest (i.e., the one created by DMForestTemplate()), then this can
454: cause a memory leak. This method is used by subtypes of DMForest when automatically constructing mesh hierarchies.
456: Logically collective on dm
458: Input Parameters:
459: + dm - the forest
460: - purpose - the adaptivity purpose
462: Level: advanced
464: .seealso: DMForestTemplate(), DMForestSetAdaptivityForest(), DMForestGetAdaptivityForest(), DMAdaptFlag
465: @*/
466: PetscErrorCode DMForestSetAdaptivityPurpose(DM dm, DMAdaptFlag purpose)
467: {
468: DM_Forest *forest;
470: forest = (DM_Forest*) dm->data;
472: if (purpose != forest->adaptPurpose) {
473: DM adapt;
475: DMForestGetAdaptivityForest(dm,&adapt);
476: PetscObjectReference((PetscObject)adapt);
477: DMForestSetAdaptivityForest(dm,NULL);
479: forest->adaptPurpose = purpose;
481: DMForestSetAdaptivityForest(dm,adapt);
482: DMDestroy(&adapt);
483: }
484: return 0;
485: }
487: /*@
488: DMForestGetAdaptivityPurpose - Get whether the current DM is being adapted from its source (set with
489: DMForestSetAdaptivityForest()) for the purpose of refinement (DM_ADAPT_REFINE), coarsening (DM_ADAPT_COARSEN),
490: coarsening only the last level (DM_ADAPT_COARSEN_LAST) or undefined (DM_ADAPT_DETERMINE).
491: This only matters for the purposes of reference counting: during DMDestroy(), cyclic
492: references can be found between DMs only if the cyclic reference is due to a fine/coarse relationship (see
493: DMSetFineDM()/DMSetCoarseDM()). If the purpose is not refinement or coarsening, and the user does not maintain a
494: reference to the post-adaptation forest (i.e., the one created by DMForestTemplate()), then this can cause a memory
495: leak. This method is used by subtypes of DMForest when automatically constructing mesh hierarchies.
497: Not collective
499: Input Parameter:
500: . dm - the forest
502: Output Parameter:
503: . purpose - the adaptivity purpose
505: Level: advanced
507: .seealso: DMForestTemplate(), DMForestSetAdaptivityForest(), DMForestGetAdaptivityForest(), DMAdaptFlag
508: @*/
509: PetscErrorCode DMForestGetAdaptivityPurpose(DM dm, DMAdaptFlag *purpose)
510: {
511: DM_Forest *forest;
513: forest = (DM_Forest*) dm->data;
514: *purpose = forest->adaptPurpose;
515: return 0;
516: }
518: /*@
519: DMForestSetAdjacencyDimension - During the pre-setup phase, set the dimension of interface points that determine
520: cell adjacency (for the purposes of partitioning and overlap).
522: Logically collective on dm
524: Input Parameters:
525: + dm - the forest
526: - adjDim - default 0 (i.e., vertices determine adjacency)
528: Level: intermediate
530: .seealso: DMForestGetAdjacencyDimension(), DMForestSetAdjacencyCodimension(), DMForestSetPartitionOverlap()
531: @*/
532: PetscErrorCode DMForestSetAdjacencyDimension(DM dm, PetscInt adjDim)
533: {
534: PetscInt dim;
535: DM_Forest *forest = (DM_Forest*) dm->data;
540: DMGetDimension(dm,&dim);
542: forest->adjDim = adjDim;
543: return 0;
544: }
546: /*@
547: DMForestSetAdjacencyCodimension - Like DMForestSetAdjacencyDimension(), but specified as a co-dimension (so that,
548: e.g., adjacency based on facets can be specified by codimension 1 in all cases)
550: Logically collective on dm
552: Input Parameters:
553: + dm - the forest
554: - adjCodim - default isthe dimension of the forest (see DMGetDimension()), since this is the codimension of vertices
556: Level: intermediate
558: .seealso: DMForestGetAdjacencyCodimension(), DMForestSetAdjacencyDimension()
559: @*/
560: PetscErrorCode DMForestSetAdjacencyCodimension(DM dm, PetscInt adjCodim)
561: {
562: PetscInt dim;
565: DMGetDimension(dm,&dim);
566: DMForestSetAdjacencyDimension(dm,dim-adjCodim);
567: return 0;
568: }
570: /*@
571: DMForestGetAdjacencyDimension - Get the dimension of interface points that determine cell adjacency (for the
572: purposes of partitioning and overlap).
574: Not collective
576: Input Parameter:
577: . dm - the forest
579: Output Parameter:
580: . adjDim - default 0 (i.e., vertices determine adjacency)
582: Level: intermediate
584: .seealso: DMForestSetAdjacencyDimension(), DMForestGetAdjacencyCodimension(), DMForestSetPartitionOverlap()
585: @*/
586: PetscErrorCode DMForestGetAdjacencyDimension(DM dm, PetscInt *adjDim)
587: {
588: DM_Forest *forest = (DM_Forest*) dm->data;
592: *adjDim = forest->adjDim;
593: return 0;
594: }
596: /*@
597: DMForestGetAdjacencyCodimension - Like DMForestGetAdjacencyDimension(), but specified as a co-dimension (so that,
598: e.g., adjacency based on facets can be specified by codimension 1 in all cases)
600: Not collective
602: Input Parameter:
603: . dm - the forest
605: Output Parameter:
606: . adjCodim - default isthe dimension of the forest (see DMGetDimension()), since this is the codimension of vertices
608: Level: intermediate
610: .seealso: DMForestSetAdjacencyCodimension(), DMForestGetAdjacencyDimension()
611: @*/
612: PetscErrorCode DMForestGetAdjacencyCodimension(DM dm, PetscInt *adjCodim)
613: {
614: DM_Forest *forest = (DM_Forest*) dm->data;
615: PetscInt dim;
619: DMGetDimension(dm,&dim);
620: *adjCodim = dim - forest->adjDim;
621: return 0;
622: }
624: /*@
625: DMForestSetPartitionOverlap - During the pre-setup phase, set the amount of cell-overlap present in parallel
626: partitions of a forest, with values > 0 indicating subdomains that are expanded by that many iterations of adding
627: adjacent cells
629: Logically collective on dm
631: Input Parameters:
632: + dm - the forest
633: - overlap - default 0
635: Level: intermediate
637: .seealso: DMForestGetPartitionOverlap(), DMForestSetAdjacencyDimension(), DMForestSetAdjacencyCodimension()
638: @*/
639: PetscErrorCode DMForestSetPartitionOverlap(DM dm, PetscInt overlap)
640: {
641: DM_Forest *forest = (DM_Forest*) dm->data;
646: forest->overlap = overlap;
647: return 0;
648: }
650: /*@
651: DMForestGetPartitionOverlap - Get the amount of cell-overlap present in parallel partitions of a forest, with values
652: > 0 indicating subdomains that are expanded by that many iterations of adding adjacent cells
654: Not collective
656: Input Parameter:
657: . dm - the forest
659: Output Parameter:
660: . overlap - default 0
662: Level: intermediate
664: .seealso: DMForestGetPartitionOverlap(), DMForestSetAdjacencyDimension(), DMForestSetAdjacencyCodimension()
665: @*/
666: PetscErrorCode DMForestGetPartitionOverlap(DM dm, PetscInt *overlap)
667: {
668: DM_Forest *forest = (DM_Forest*) dm->data;
672: *overlap = forest->overlap;
673: return 0;
674: }
676: /*@
677: DMForestSetMinimumRefinement - During the pre-setup phase, set the minimum level of refinement (relative to the base
678: DM, see DMForestGetBaseDM()) allowed in the forest. If the forest is being created by coarsening a previous forest
679: (see DMForestGetAdaptivityForest()) this limits the amount of coarsening.
681: Logically collective on dm
683: Input Parameters:
684: + dm - the forest
685: - minRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
687: Level: intermediate
689: .seealso: DMForestGetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestSetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityForest()
690: @*/
691: PetscErrorCode DMForestSetMinimumRefinement(DM dm, PetscInt minRefinement)
692: {
693: DM_Forest *forest = (DM_Forest*) dm->data;
697: forest->minRefinement = minRefinement;
698: return 0;
699: }
701: /*@
702: DMForestGetMinimumRefinement - Get the minimum level of refinement (relative to the base DM, see
703: DMForestGetBaseDM()) allowed in the forest. If the forest is being created by coarsening a previous forest (see
704: DMForestGetAdaptivityForest()), this limits the amount of coarsening.
706: Not collective
708: Input Parameter:
709: . dm - the forest
711: Output Parameter:
712: . minRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
714: Level: intermediate
716: .seealso: DMForestSetMinimumRefinement(), DMForestGetMaximumRefinement(), DMForestGetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityForest()
717: @*/
718: PetscErrorCode DMForestGetMinimumRefinement(DM dm, PetscInt *minRefinement)
719: {
720: DM_Forest *forest = (DM_Forest*) dm->data;
724: *minRefinement = forest->minRefinement;
725: return 0;
726: }
728: /*@
729: DMForestSetInitialRefinement - During the pre-setup phase, set the initial level of refinement (relative to the base
730: DM, see DMForestGetBaseDM()) allowed in the forest.
732: Logically collective on dm
734: Input Parameters:
735: + dm - the forest
736: - initefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
738: Level: intermediate
740: .seealso: DMForestSetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestGetBaseDM()
741: @*/
742: PetscErrorCode DMForestSetInitialRefinement(DM dm, PetscInt initRefinement)
743: {
744: DM_Forest *forest = (DM_Forest*) dm->data;
748: forest->initRefinement = initRefinement;
749: return 0;
750: }
752: /*@
753: DMForestGetInitialRefinement - Get the initial level of refinement (relative to the base DM, see
754: DMForestGetBaseDM()) allowed in the forest.
756: Not collective
758: Input Parameter:
759: . dm - the forest
761: Output Parameter:
762: . initRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
764: Level: intermediate
766: .seealso: DMForestSetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestGetBaseDM()
767: @*/
768: PetscErrorCode DMForestGetInitialRefinement(DM dm, PetscInt *initRefinement)
769: {
770: DM_Forest *forest = (DM_Forest*) dm->data;
774: *initRefinement = forest->initRefinement;
775: return 0;
776: }
778: /*@
779: DMForestSetMaximumRefinement - During the pre-setup phase, set the maximum level of refinement (relative to the base
780: DM, see DMForestGetBaseDM()) allowed in the forest. If the forest is being created by refining a previous forest
781: (see DMForestGetAdaptivityForest()), this limits the amount of refinement.
783: Logically collective on dm
785: Input Parameters:
786: + dm - the forest
787: - maxRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
789: Level: intermediate
791: .seealso: DMForestGetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestSetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityDM()
792: @*/
793: PetscErrorCode DMForestSetMaximumRefinement(DM dm, PetscInt maxRefinement)
794: {
795: DM_Forest *forest = (DM_Forest*) dm->data;
799: forest->maxRefinement = maxRefinement;
800: return 0;
801: }
803: /*@
804: DMForestGetMaximumRefinement - Get the maximum level of refinement (relative to the base DM, see
805: DMForestGetBaseDM()) allowed in the forest. If the forest is being created by refining a previous forest (see
806: DMForestGetAdaptivityForest()), this limits the amount of refinement.
808: Not collective
810: Input Parameter:
811: . dm - the forest
813: Output Parameter:
814: . maxRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
816: Level: intermediate
818: .seealso: DMForestSetMaximumRefinement(), DMForestGetMinimumRefinement(), DMForestGetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityForest()
819: @*/
820: PetscErrorCode DMForestGetMaximumRefinement(DM dm, PetscInt *maxRefinement)
821: {
822: DM_Forest *forest = (DM_Forest*) dm->data;
826: *maxRefinement = forest->maxRefinement;
827: return 0;
828: }
830: /*@C
831: DMForestSetAdaptivityStrategy - During the pre-setup phase, set the strategy for combining adaptivity labels from multiple processes.
832: Subtypes of DMForest may define their own strategies. Two default strategies are DMFORESTADAPTALL, which indicates that all processes must agree
833: for a refinement/coarsening flag to be valid, and DMFORESTADAPTANY, which indicates that only one process needs to
834: specify refinement/coarsening.
836: Logically collective on dm
838: Input Parameters:
839: + dm - the forest
840: - adaptStrategy - default DMFORESTADAPTALL
842: Level: advanced
844: .seealso: DMForestGetAdaptivityStrategy()
845: @*/
846: PetscErrorCode DMForestSetAdaptivityStrategy(DM dm, DMForestAdaptivityStrategy adaptStrategy)
847: {
848: DM_Forest *forest = (DM_Forest*) dm->data;
851: PetscFree(forest->adaptStrategy);
852: PetscStrallocpy((const char*) adaptStrategy,(char**)&forest->adaptStrategy);
853: return 0;
854: }
856: /*@C
857: DMForestSetAdaptivityStrategy - Get the strategy for combining adaptivity labels from multiple processes. Subtypes
858: of DMForest may define their own strategies. Two default strategies are DMFORESTADAPTALL, which indicates that all
859: processes must agree for a refinement/coarsening flag to be valid, and DMFORESTADAPTANY, which indicates that only
860: one process needs to specify refinement/coarsening.
862: Not collective
864: Input Parameter:
865: . dm - the forest
867: Output Parameter:
868: . adaptStrategy - the adaptivity strategy (default DMFORESTADAPTALL)
870: Level: advanced
872: .seealso: DMForestSetAdaptivityStrategy()
873: @*/
874: PetscErrorCode DMForestGetAdaptivityStrategy(DM dm, DMForestAdaptivityStrategy *adaptStrategy)
875: {
876: DM_Forest *forest = (DM_Forest*) dm->data;
880: *adaptStrategy = forest->adaptStrategy;
881: return 0;
882: }
884: /*@
885: DMForestGetAdaptivitySuccess - Return whether the requested adaptation (refinement, coarsening, repartitioning,
886: etc.) was successful. PETSC_FALSE indicates that the post-adaptation forest is the same as the pre-adpatation
887: forest. A requested adaptation may have been unsuccessful if, for example, the requested refinement would have
888: exceeded the maximum refinement level.
890: Collective on dm
892: Input Parameter:
894: . dm - the post-adaptation forest
896: Output Parameter:
898: . success - PETSC_TRUE if the post-adaptation forest is different from the pre-adaptation forest.
900: Level: intermediate
902: .see
903: @*/
904: PetscErrorCode DMForestGetAdaptivitySuccess(DM dm, PetscBool *success)
905: {
906: DM_Forest *forest;
910: forest = (DM_Forest *) dm->data;
911: (forest->getadaptivitysuccess)(dm,success);
912: return 0;
913: }
915: /*@
916: DMForestSetComputeAdaptivitySF - During the pre-setup phase, set whether transfer PetscSFs should be computed
917: 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().
919: Logically collective on dm
921: Input Parameters:
922: + dm - the post-adaptation forest
923: - computeSF - default PETSC_TRUE
925: Level: advanced
927: .seealso: DMForestGetComputeAdaptivitySF(), DMForestGetAdaptivitySF()
928: @*/
929: PetscErrorCode DMForestSetComputeAdaptivitySF(DM dm, PetscBool computeSF)
930: {
931: DM_Forest *forest;
935: forest = (DM_Forest*) dm->data;
936: forest->computeAdaptSF = computeSF;
937: return 0;
938: }
940: PetscErrorCode DMForestTransferVec(DM dmIn, Vec vecIn, DM dmOut, Vec vecOut, PetscBool useBCs, PetscReal time)
941: {
942: DM_Forest *forest;
948: forest = (DM_Forest *) dmIn->data;
950: (forest->transfervec)(dmIn,vecIn,dmOut,vecOut,useBCs,time);
951: return 0;
952: }
954: PetscErrorCode DMForestTransferVecFromBase(DM dm, Vec vecIn, Vec vecOut)
955: {
956: DM_Forest *forest;
961: forest = (DM_Forest *) dm->data;
963: (forest->transfervecfrombase)(dm,vecIn,vecOut);
964: return 0;
965: }
967: /*@
968: DMForestGetComputeAdaptivitySF - Get whether transfer PetscSFs should be computed relating the cells of the
969: pre-adaptation forest to the post-adaptiation forest. After DMSetUp() is called, these transfer PetscSFs can be
970: accessed with DMForestGetAdaptivitySF().
972: Not collective
974: Input Parameter:
975: . dm - the post-adaptation forest
977: Output Parameter:
978: . computeSF - default PETSC_TRUE
980: Level: advanced
982: .seealso: DMForestSetComputeAdaptivitySF(), DMForestGetAdaptivitySF()
983: @*/
984: PetscErrorCode DMForestGetComputeAdaptivitySF(DM dm, PetscBool *computeSF)
985: {
986: DM_Forest *forest;
989: forest = (DM_Forest*) dm->data;
990: *computeSF = forest->computeAdaptSF;
991: return 0;
992: }
994: /*@
995: DMForestGetAdaptivitySF - Get PetscSFs that relate the pre-adaptation forest to the post-adaptation forest.
996: Adaptation can be any combination of refinement, coarsening, repartition, and change of overlap, so there may be
997: some cells of the pre-adaptation that are parents of post-adaptation cells, and vice versa. Therefore there are two
998: PetscSFs: one that relates pre-adaptation coarse cells to post-adaptation fine cells, and one that relates
999: pre-adaptation fine cells to post-adaptation coarse cells.
1001: Not collective
1003: Input Parameter:
1004: dm - the post-adaptation forest
1006: Output Parameter:
1007: preCoarseToFine - pre-adaptation coarse cells to post-adaptation fine cells: BCast goes from pre- to post-
1008: coarseToPreFine - post-adaptation coarse cells to pre-adaptation fine cells: BCast goes from post- to pre-
1010: Level: advanced
1012: .seealso: DMForestGetComputeAdaptivitySF(), DMForestSetComputeAdaptivitySF()
1013: @*/
1014: PetscErrorCode DMForestGetAdaptivitySF(DM dm, PetscSF *preCoarseToFine, PetscSF *coarseToPreFine)
1015: {
1016: DM_Forest *forest;
1019: DMSetUp(dm);
1020: forest = (DM_Forest*) dm->data;
1021: if (preCoarseToFine) *preCoarseToFine = forest->preCoarseToFine;
1022: if (coarseToPreFine) *coarseToPreFine = forest->coarseToPreFine;
1023: return 0;
1024: }
1026: /*@
1027: DMForestSetGradeFactor - During the pre-setup phase, set the desired amount of grading in the mesh, e.g. give 2 to
1028: indicate that the diameter of neighboring cells should differ by at most a factor of 2. Subtypes of DMForest may
1029: only support one particular choice of grading factor.
1031: Logically collective on dm
1033: Input Parameters:
1034: + dm - the forest
1035: - grade - the grading factor
1037: Level: advanced
1039: .seealso: DMForestGetGradeFactor()
1040: @*/
1041: PetscErrorCode DMForestSetGradeFactor(DM dm, PetscInt grade)
1042: {
1043: DM_Forest *forest = (DM_Forest*) dm->data;
1047: forest->gradeFactor = grade;
1048: return 0;
1049: }
1051: /*@
1052: DMForestGetGradeFactor - Get the desired amount of grading in the mesh, e.g. give 2 to indicate that the diameter of
1053: neighboring cells should differ by at most a factor of 2. Subtypes of DMForest may only support one particular
1054: choice of grading factor.
1056: Not collective
1058: Input Parameter:
1059: . dm - the forest
1061: Output Parameter:
1062: . grade - the grading factor
1064: Level: advanced
1066: .seealso: DMForestSetGradeFactor()
1067: @*/
1068: PetscErrorCode DMForestGetGradeFactor(DM dm, PetscInt *grade)
1069: {
1070: DM_Forest *forest = (DM_Forest*) dm->data;
1074: *grade = forest->gradeFactor;
1075: return 0;
1076: }
1078: /*@
1079: DMForestSetCellWeightFactor - During the pre-setup phase, set the factor by which the level of refinement changes
1080: the cell weight (see DMForestSetCellWeights()) when calculating partitions. The final weight of a cell will be
1081: (cellWeight) * (weightFactor^refinementLevel). A factor of 1 indicates that the weight of a cell does not depend on
1082: its level; a factor of 2, for example, might be appropriate for sub-cycling time-stepping methods, when the
1083: computation associated with a cell is multiplied by a factor of 2 for each additional level of refinement.
1085: Logically collective on dm
1087: Input Parameters:
1088: + dm - the forest
1089: - weightsFactors - default 1.
1091: Level: advanced
1093: .seealso: DMForestGetCellWeightFactor(), DMForestSetCellWeights()
1094: @*/
1095: PetscErrorCode DMForestSetCellWeightFactor(DM dm, PetscReal weightsFactor)
1096: {
1097: DM_Forest *forest = (DM_Forest*) dm->data;
1101: forest->weightsFactor = weightsFactor;
1102: return 0;
1103: }
1105: /*@
1106: DMForestGetCellWeightFactor - Get the factor by which the level of refinement changes the cell weight (see
1107: DMForestSetCellWeights()) when calculating partitions. The final weight of a cell will be (cellWeight) *
1108: (weightFactor^refinementLevel). A factor of 1 indicates that the weight of a cell does not depend on its level; a
1109: factor of 2, for example, might be appropriate for sub-cycling time-stepping methods, when the computation
1110: associated with a cell is multiplied by a factor of 2 for each additional level of refinement.
1112: Not collective
1114: Input Parameter:
1115: . dm - the forest
1117: Output Parameter:
1118: . weightsFactors - default 1.
1120: Level: advanced
1122: .seealso: DMForestSetCellWeightFactor(), DMForestSetCellWeights()
1123: @*/
1124: PetscErrorCode DMForestGetCellWeightFactor(DM dm, PetscReal *weightsFactor)
1125: {
1126: DM_Forest *forest = (DM_Forest*) dm->data;
1130: *weightsFactor = forest->weightsFactor;
1131: return 0;
1132: }
1134: /*@
1135: DMForestGetCellChart - After the setup phase, get the local half-open interval of the chart of cells on this process
1137: Not collective
1139: Input Parameter:
1140: . dm - the forest
1142: Output Parameters:
1143: + cStart - the first cell on this process
1144: - cEnd - one after the final cell on this process
1146: Level: intermediate
1148: .seealso: DMForestGetCellSF()
1149: @*/
1150: PetscErrorCode DMForestGetCellChart(DM dm, PetscInt *cStart, PetscInt *cEnd)
1151: {
1152: DM_Forest *forest = (DM_Forest*) dm->data;
1157: if (((forest->cStart == PETSC_DETERMINE) || (forest->cEnd == PETSC_DETERMINE)) && forest->createcellchart) {
1158: forest->createcellchart(dm,&forest->cStart,&forest->cEnd);
1159: }
1160: *cStart = forest->cStart;
1161: *cEnd = forest->cEnd;
1162: return 0;
1163: }
1165: /*@
1166: DMForestGetCellSF - After the setup phase, get the PetscSF for overlapping cells between processes
1168: Not collective
1170: Input Parameter:
1171: . dm - the forest
1173: Output Parameter:
1174: . cellSF - the PetscSF
1176: Level: intermediate
1178: .seealso: DMForestGetCellChart()
1179: @*/
1180: PetscErrorCode DMForestGetCellSF(DM dm, PetscSF *cellSF)
1181: {
1182: DM_Forest *forest = (DM_Forest*) dm->data;
1186: if ((!forest->cellSF) && forest->createcellsf) {
1187: forest->createcellsf(dm,&forest->cellSF);
1188: }
1189: *cellSF = forest->cellSF;
1190: return 0;
1191: }
1193: /*@C
1194: DMForestSetAdaptivityLabel - During the pre-setup phase, set the label of the pre-adaptation forest (see
1195: DMForestGetAdaptivityForest()) that holds the adaptation flags (refinement, coarsening, or some combination). The
1196: interpretation of the label values is up to the subtype of DMForest, but DM_ADAPT_DETERMINE, DM_ADAPT_KEEP,
1197: DM_ADAPT_REFINE, and DM_ADAPT_COARSEN have been reserved as choices that should be accepted by all subtypes.
1199: Logically collective on dm
1201: Input Parameters:
1202: - dm - the forest
1203: + adaptLabel - the label in the pre-adaptation forest
1205: Level: intermediate
1207: .seealso DMForestGetAdaptivityLabel()
1208: @*/
1209: PetscErrorCode DMForestSetAdaptivityLabel(DM dm, DMLabel adaptLabel)
1210: {
1211: DM_Forest *forest = (DM_Forest*) dm->data;
1215: PetscObjectReference((PetscObject)adaptLabel);
1216: DMLabelDestroy(&forest->adaptLabel);
1217: forest->adaptLabel = adaptLabel;
1218: return 0;
1219: }
1221: /*@C
1222: DMForestGetAdaptivityLabel - Get the label of the pre-adaptation forest (see DMForestGetAdaptivityForest()) that
1223: holds the adaptation flags (refinement, coarsening, or some combination). The interpretation of the label values is
1224: up to the subtype of DMForest, but DM_ADAPT_DETERMINE, DM_ADAPT_KEEP, DM_ADAPT_REFINE, and DM_ADAPT_COARSEN have
1225: been reserved as choices that should be accepted by all subtypes.
1227: Not collective
1229: Input Parameter:
1230: . dm - the forest
1232: Output Parameter:
1233: . adaptLabel - the name of the label in the pre-adaptation forest
1235: Level: intermediate
1237: .seealso DMForestSetAdaptivityLabel()
1238: @*/
1239: PetscErrorCode DMForestGetAdaptivityLabel(DM dm, DMLabel *adaptLabel)
1240: {
1241: DM_Forest *forest = (DM_Forest*) dm->data;
1244: *adaptLabel = forest->adaptLabel;
1245: return 0;
1246: }
1248: /*@
1249: DMForestSetCellWeights - Set the weights assigned to each of the cells (see DMForestGetCellChart()) of the current
1250: process: weights are used to determine parallel partitioning. Partitions will be created so that each process's
1251: ratio of weight to capacity (see DMForestSetWeightCapacity()) is roughly equal. If NULL, each cell receives a weight
1252: of 1.
1254: Logically collective on dm
1256: Input Parameters:
1257: + dm - the forest
1258: . weights - the array of weights for all cells, or NULL to indicate each cell has weight 1.
1259: - copyMode - how weights should reference weights
1261: Level: advanced
1263: .seealso: DMForestGetCellWeights(), DMForestSetWeightCapacity()
1264: @*/
1265: PetscErrorCode DMForestSetCellWeights(DM dm, PetscReal weights[], PetscCopyMode copyMode)
1266: {
1267: DM_Forest *forest = (DM_Forest*) dm->data;
1268: PetscInt cStart, cEnd;
1271: DMForestGetCellChart(dm,&cStart,&cEnd);
1273: if (copyMode == PETSC_COPY_VALUES) {
1274: if (forest->cellWeightsCopyMode != PETSC_OWN_POINTER || forest->cellWeights == weights) {
1275: PetscMalloc1(cEnd-cStart,&forest->cellWeights);
1276: }
1277: PetscArraycpy(forest->cellWeights,weights,cEnd-cStart);
1278: forest->cellWeightsCopyMode = PETSC_OWN_POINTER;
1279: return 0;
1280: }
1281: if (forest->cellWeightsCopyMode == PETSC_OWN_POINTER) {
1282: PetscFree(forest->cellWeights);
1283: }
1284: forest->cellWeights = weights;
1285: forest->cellWeightsCopyMode = copyMode;
1286: return 0;
1287: }
1289: /*@
1290: DMForestGetCellWeights - Get the weights assigned to each of the cells (see DMForestGetCellChart()) of the current
1291: process: weights are used to determine parallel partitioning. Partitions will be created so that each process's
1292: ratio of weight to capacity (see DMForestSetWeightCapacity()) is roughly equal. If NULL, each cell receives a weight
1293: of 1.
1295: Not collective
1297: Input Parameter:
1298: . dm - the forest
1300: Output Parameter:
1301: . weights - the array of weights for all cells, or NULL to indicate each cell has weight 1.
1303: Level: advanced
1305: .seealso: DMForestSetCellWeights(), DMForestSetWeightCapacity()
1306: @*/
1307: PetscErrorCode DMForestGetCellWeights(DM dm, PetscReal **weights)
1308: {
1309: DM_Forest *forest = (DM_Forest*) dm->data;
1313: *weights = forest->cellWeights;
1314: return 0;
1315: }
1317: /*@
1318: DMForestSetWeightCapacity - During the pre-setup phase, set the capacity of the current process when repartitioning
1319: a pre-adaptation forest (see DMForestGetAdaptivityForest()). After partitioning, the ratio of the weight of each
1320: process's cells to the process's capacity will be roughly equal for all processes. A capacity of 0 indicates that
1321: the current process should not have any cells after repartitioning.
1323: Logically Collective on dm
1325: Input parameters:
1326: + dm - the forest
1327: - capacity - this process's capacity
1329: Level: advanced
1331: .seealso DMForestGetWeightCapacity(), DMForestSetCellWeights(), DMForestSetCellWeightFactor()
1332: @*/
1333: PetscErrorCode DMForestSetWeightCapacity(DM dm, PetscReal capacity)
1334: {
1335: DM_Forest *forest = (DM_Forest*) dm->data;
1340: forest->weightCapacity = capacity;
1341: return 0;
1342: }
1344: /*@
1345: DMForestGetWeightCapacity - Set the capacity of the current process when repartitioning a pre-adaptation forest (see
1346: DMForestGetAdaptivityForest()). After partitioning, the ratio of the weight of each process's cells to the
1347: process's capacity will be roughly equal for all processes. A capacity of 0 indicates that the current process
1348: should not have any cells after repartitioning.
1350: Not collective
1352: Input parameter:
1353: . dm - the forest
1355: Output parameter:
1356: . capacity - this process's capacity
1358: Level: advanced
1360: .seealso DMForestSetWeightCapacity(), DMForestSetCellWeights(), DMForestSetCellWeightFactor()
1361: @*/
1362: PetscErrorCode DMForestGetWeightCapacity(DM dm, PetscReal *capacity)
1363: {
1364: DM_Forest *forest = (DM_Forest*) dm->data;
1368: *capacity = forest->weightCapacity;
1369: return 0;
1370: }
1372: PETSC_EXTERN PetscErrorCode DMSetFromOptions_Forest(PetscOptionItems *PetscOptionsObject,DM dm)
1373: {
1374: PetscBool flg, flg1, flg2, flg3, flg4;
1375: DMForestTopology oldTopo;
1376: char stringBuffer[256];
1377: PetscViewer viewer;
1378: PetscViewerFormat format;
1379: PetscInt adjDim, adjCodim, overlap, minRefinement, initRefinement, maxRefinement, grade;
1380: PetscReal weightsFactor;
1381: DMForestAdaptivityStrategy adaptStrategy;
1384: DMForestGetTopology(dm, &oldTopo);
1385: PetscOptionsHead(PetscOptionsObject,"DMForest Options");
1386: PetscOptionsString("-dm_forest_topology","the topology of the forest's base mesh","DMForestSetTopology",oldTopo,stringBuffer,sizeof(stringBuffer),&flg1);
1387: PetscOptionsViewer("-dm_forest_base_dm","load the base DM from a viewer specification","DMForestSetBaseDM",&viewer,&format,&flg2);
1388: PetscOptionsViewer("-dm_forest_coarse_forest","load the coarse forest from a viewer specification","DMForestSetCoarseForest",&viewer,&format,&flg3);
1389: PetscOptionsViewer("-dm_forest_fine_forest","load the fine forest from a viewer specification","DMForestSetFineForest",&viewer,&format,&flg4);
1391: if (flg1) {
1392: DMForestSetTopology(dm,(DMForestTopology)stringBuffer);
1393: DMForestSetBaseDM(dm,NULL);
1394: DMForestSetAdaptivityForest(dm,NULL);
1395: }
1396: if (flg2) {
1397: DM base;
1399: DMCreate(PetscObjectComm((PetscObject)dm),&base);
1400: PetscViewerPushFormat(viewer,format);
1401: DMLoad(base,viewer);
1402: PetscViewerDestroy(&viewer);
1403: DMForestSetBaseDM(dm,base);
1404: DMDestroy(&base);
1405: DMForestSetTopology(dm,NULL);
1406: DMForestSetAdaptivityForest(dm,NULL);
1407: }
1408: if (flg3) {
1409: DM coarse;
1411: DMCreate(PetscObjectComm((PetscObject)dm),&coarse);
1412: PetscViewerPushFormat(viewer,format);
1413: DMLoad(coarse,viewer);
1414: PetscViewerDestroy(&viewer);
1415: DMForestSetAdaptivityForest(dm,coarse);
1416: DMDestroy(&coarse);
1417: DMForestSetTopology(dm,NULL);
1418: DMForestSetBaseDM(dm,NULL);
1419: }
1420: if (flg4) {
1421: DM fine;
1423: DMCreate(PetscObjectComm((PetscObject)dm),&fine);
1424: PetscViewerPushFormat(viewer,format);
1425: DMLoad(fine,viewer);
1426: PetscViewerDestroy(&viewer);
1427: DMForestSetAdaptivityForest(dm,fine);
1428: DMDestroy(&fine);
1429: DMForestSetTopology(dm,NULL);
1430: DMForestSetBaseDM(dm,NULL);
1431: }
1432: DMForestGetAdjacencyDimension(dm,&adjDim);
1433: PetscOptionsBoundedInt("-dm_forest_adjacency_dimension","set the dimension of points that define adjacency in the forest","DMForestSetAdjacencyDimension",adjDim,&adjDim,&flg,0);
1434: if (flg) {
1435: DMForestSetAdjacencyDimension(dm,adjDim);
1436: } else {
1437: DMForestGetAdjacencyCodimension(dm,&adjCodim);
1438: PetscOptionsBoundedInt("-dm_forest_adjacency_codimension","set the codimension of points that define adjacency in the forest","DMForestSetAdjacencyCodimension",adjCodim,&adjCodim,&flg,1);
1439: if (flg) {
1440: DMForestSetAdjacencyCodimension(dm,adjCodim);
1441: }
1442: }
1443: DMForestGetPartitionOverlap(dm,&overlap);
1444: PetscOptionsBoundedInt("-dm_forest_partition_overlap","set the degree of partition overlap","DMForestSetPartitionOverlap",overlap,&overlap,&flg,0);
1445: if (flg) {
1446: DMForestSetPartitionOverlap(dm,overlap);
1447: }
1448: #if 0
1449: 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);
1450: if (flg) {
1451: DMForestSetMinimumRefinement(dm,minRefinement);
1452: DMForestSetInitialRefinement(dm,minRefinement);
1453: }
1454: PetscOptionsBoundedInt("-dm_refine_hierarchy","equivalent to -dm_forest_set_minimum_refinement 0 and -dm_forest_set_initial_refinement",NULL,initRefinement,&initRefinement,&flg,0);
1455: if (flg) {
1456: DMForestSetMinimumRefinement(dm,0);
1457: DMForestSetInitialRefinement(dm,initRefinement);
1458: }
1459: #endif
1460: DMForestGetMinimumRefinement(dm,&minRefinement);
1461: PetscOptionsBoundedInt("-dm_forest_minimum_refinement","set the minimum level of refinement in the forest","DMForestSetMinimumRefinement",minRefinement,&minRefinement,&flg,0);
1462: if (flg) {
1463: DMForestSetMinimumRefinement(dm,minRefinement);
1464: }
1465: DMForestGetInitialRefinement(dm,&initRefinement);
1466: PetscOptionsBoundedInt("-dm_forest_initial_refinement","set the initial level of refinement in the forest","DMForestSetInitialRefinement",initRefinement,&initRefinement,&flg,0);
1467: if (flg) {
1468: DMForestSetInitialRefinement(dm,initRefinement);
1469: }
1470: DMForestGetMaximumRefinement(dm,&maxRefinement);
1471: PetscOptionsBoundedInt("-dm_forest_maximum_refinement","set the maximum level of refinement in the forest","DMForestSetMaximumRefinement",maxRefinement,&maxRefinement,&flg,0);
1472: if (flg) {
1473: DMForestSetMaximumRefinement(dm,maxRefinement);
1474: }
1475: DMForestGetAdaptivityStrategy(dm,&adaptStrategy);
1476: PetscOptionsString("-dm_forest_adaptivity_strategy","the forest's adaptivity-flag resolution strategy","DMForestSetAdaptivityStrategy",adaptStrategy,stringBuffer,sizeof(stringBuffer),&flg);
1477: if (flg) {
1478: DMForestSetAdaptivityStrategy(dm,(DMForestAdaptivityStrategy)stringBuffer);
1479: }
1480: DMForestGetGradeFactor(dm,&grade);
1481: PetscOptionsBoundedInt("-dm_forest_grade_factor","grade factor between neighboring cells","DMForestSetGradeFactor",grade,&grade,&flg,0);
1482: if (flg) {
1483: DMForestSetGradeFactor(dm,grade);
1484: }
1485: DMForestGetCellWeightFactor(dm,&weightsFactor);
1486: PetscOptionsReal("-dm_forest_cell_weight_factor","multiplying weight factor for cell refinement","DMForestSetCellWeightFactor",weightsFactor,&weightsFactor,&flg);
1487: if (flg) {
1488: DMForestSetCellWeightFactor(dm,weightsFactor);
1489: }
1490: PetscOptionsTail();
1491: return 0;
1492: }
1494: PetscErrorCode DMCreateSubDM_Forest(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1495: {
1496: if (subdm) DMClone(dm, subdm);
1497: DMCreateSectionSubDM(dm, numFields, fields, is, subdm);
1498: return 0;
1499: }
1501: PetscErrorCode DMRefine_Forest(DM dm, MPI_Comm comm, DM *dmRefined)
1502: {
1503: DMLabel refine;
1504: DM fineDM;
1506: DMGetFineDM(dm,&fineDM);
1507: if (fineDM) {
1508: PetscObjectReference((PetscObject)fineDM);
1509: *dmRefined = fineDM;
1510: return 0;
1511: }
1512: DMForestTemplate(dm,comm,dmRefined);
1513: DMGetLabel(dm,"refine",&refine);
1514: if (!refine) {
1515: DMLabelCreate(PETSC_COMM_SELF, "refine",&refine);
1516: DMLabelSetDefaultValue(refine,DM_ADAPT_REFINE);
1517: } else {
1518: PetscObjectReference((PetscObject) refine);
1519: }
1520: DMForestSetAdaptivityLabel(*dmRefined,refine);
1521: DMLabelDestroy(&refine);
1522: return 0;
1523: }
1525: PetscErrorCode DMCoarsen_Forest(DM dm, MPI_Comm comm, DM *dmCoarsened)
1526: {
1527: DMLabel coarsen;
1528: DM coarseDM;
1530: {
1531: PetscMPIInt mpiComparison;
1532: MPI_Comm dmcomm = PetscObjectComm((PetscObject)dm);
1534: MPI_Comm_compare(comm, dmcomm, &mpiComparison);
1536: }
1537: DMGetCoarseDM(dm,&coarseDM);
1538: if (coarseDM) {
1539: PetscObjectReference((PetscObject)coarseDM);
1540: *dmCoarsened = coarseDM;
1541: return 0;
1542: }
1543: DMForestTemplate(dm,comm,dmCoarsened);
1544: DMForestSetAdaptivityPurpose(*dmCoarsened,DM_ADAPT_COARSEN);
1545: DMGetLabel(dm,"coarsen",&coarsen);
1546: if (!coarsen) {
1547: DMLabelCreate(PETSC_COMM_SELF, "coarsen",&coarsen);
1548: DMLabelSetDefaultValue(coarsen,DM_ADAPT_COARSEN);
1549: } else {
1550: PetscObjectReference((PetscObject) coarsen);
1551: }
1552: DMForestSetAdaptivityLabel(*dmCoarsened,coarsen);
1553: DMLabelDestroy(&coarsen);
1554: return 0;
1555: }
1557: PetscErrorCode DMAdaptLabel_Forest(DM dm, PETSC_UNUSED Vec metric, DMLabel label, PETSC_UNUSED DMLabel rgLabel, DM *adaptedDM)
1558: {
1559: PetscBool success;
1561: DMForestTemplate(dm,PetscObjectComm((PetscObject)dm),adaptedDM);
1562: DMForestSetAdaptivityLabel(*adaptedDM,label);
1563: DMSetUp(*adaptedDM);
1564: DMForestGetAdaptivitySuccess(*adaptedDM,&success);
1565: if (!success) {
1566: DMDestroy(adaptedDM);
1567: *adaptedDM = NULL;
1568: }
1569: return 0;
1570: }
1572: static PetscErrorCode DMInitialize_Forest(DM dm)
1573: {
1574: PetscMemzero(dm->ops,sizeof(*(dm->ops)));
1576: dm->ops->clone = DMClone_Forest;
1577: dm->ops->setfromoptions = DMSetFromOptions_Forest;
1578: dm->ops->destroy = DMDestroy_Forest;
1579: dm->ops->createsubdm = DMCreateSubDM_Forest;
1580: dm->ops->refine = DMRefine_Forest;
1581: dm->ops->coarsen = DMCoarsen_Forest;
1582: return 0;
1583: }
1585: /*MC
1587: DMFOREST = "forest" - A DM object that encapsulates a hierarchically refined mesh. Forests usually have a base DM
1588: (see DMForestGetBaseDM()), from which it is refined. The refinement and partitioning of forests is considered
1589: immutable after DMSetUp() is called. To adapt a mesh, one should call DMForestTemplate() to create a new mesh that
1590: will default to being identical to it, specify how that mesh should differ, and then calling DMSetUp() on the new
1591: mesh.
1593: To specify that a mesh should be refined or coarsened from the previous mesh, a label should be defined on the
1594: previous mesh whose values indicate which cells should be refined (DM_ADAPT_REFINE) or coarsened (DM_ADAPT_COARSEN)
1595: and how (subtypes are free to allow additional values for things like anisotropic refinement). The label should be
1596: given to the *new* mesh with DMForestSetAdaptivityLabel().
1598: Level: advanced
1600: .seealso: DMType, DMCreate(), DMSetType(), DMForestGetBaseDM(), DMForestSetBaseDM(), DMForestTemplate(), DMForestSetAdaptivityLabel()
1601: M*/
1603: PETSC_EXTERN PetscErrorCode DMCreate_Forest(DM dm)
1604: {
1605: DM_Forest *forest;
1608: PetscNewLog(dm,&forest);
1609: dm->dim = 0;
1610: dm->data = forest;
1611: forest->refct = 1;
1612: forest->data = NULL;
1613: forest->topology = NULL;
1614: forest->adapt = NULL;
1615: forest->base = NULL;
1616: forest->adaptPurpose = DM_ADAPT_DETERMINE;
1617: forest->adjDim = PETSC_DEFAULT;
1618: forest->overlap = PETSC_DEFAULT;
1619: forest->minRefinement = PETSC_DEFAULT;
1620: forest->maxRefinement = PETSC_DEFAULT;
1621: forest->initRefinement = PETSC_DEFAULT;
1622: forest->cStart = PETSC_DETERMINE;
1623: forest->cEnd = PETSC_DETERMINE;
1624: forest->cellSF = NULL;
1625: forest->adaptLabel = NULL;
1626: forest->gradeFactor = 2;
1627: forest->cellWeights = NULL;
1628: forest->cellWeightsCopyMode = PETSC_USE_POINTER;
1629: forest->weightsFactor = 1.;
1630: forest->weightCapacity = 1.;
1631: DMForestSetAdaptivityStrategy(dm,DMFORESTADAPTALL);
1632: DMInitialize_Forest(dm);
1633: return 0;
1634: }