Actual source code: forest.c
petsc-3.14.6 2021-03-30
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;
21: PetscErrorCode ierr;
24: while (link) {
25: oldLink = link;
26: PetscFree(oldLink->name);
27: link = oldLink->next;
28: PetscFree(oldLink);
29: }
30: return(0);
31: }
33: static PetscErrorCode DMForestPackageInitialize(void)
34: {
38: if (DMForestPackageInitialized) return(0);
39: DMForestPackageInitialized = PETSC_TRUE;
41: DMForestRegisterType(DMFOREST);
42: PetscRegisterFinalize(DMForestPackageFinalize);
43: return(0);
44: }
46: /*@C
47: DMForestRegisterType - Registers a DMType as a subtype of DMFOREST (so that DMIsForest() will be correct)
49: Not Collective
51: Input parameter:
52: . name - the name of the type
54: Level: advanced
56: .seealso: DMFOREST, DMIsForest()
57: @*/
58: PetscErrorCode DMForestRegisterType(DMType name)
59: {
60: DMForestTypeLink link;
61: PetscErrorCode ierr;
64: DMForestPackageInitialize();
65: PetscNew(&link);
66: PetscStrallocpy(name,&link->name);
67: link->next = DMForestTypeList;
68: DMForestTypeList = link;
69: return(0);
70: }
72: /*@
73: DMIsForest - Check whether a DM uses the DMFOREST interface for hierarchically-refined meshes
75: Not Collective
77: Input parameter:
78: . dm - the DM object
80: Output parameter:
81: . isForest - whether dm is a subtype of DMFOREST
83: Level: intermediate
85: .seealso: DMFOREST, DMForestRegisterType()
86: @*/
87: PetscErrorCode DMIsForest(DM dm, PetscBool *isForest)
88: {
89: DMForestTypeLink link = DMForestTypeList;
90: PetscErrorCode ierr;
93: while (link) {
94: PetscBool sameType;
95: PetscObjectTypeCompare((PetscObject)dm,link->name,&sameType);
96: if (sameType) {
97: *isForest = PETSC_TRUE;
98: return(0);
99: }
100: link = link->next;
101: }
102: *isForest = PETSC_FALSE;
103: return(0);
104: }
106: /*@
107: DMForestTemplate - Create a new DM that will be adapted from a source DM. The new DM reproduces the configuration
108: of the source, but is not yet setup, so that the user can then define only the ways that the new DM should differ
109: (by, e.g., refinement or repartitioning). The source DM is also set as the adaptivity source DM of the new DM (see
110: DMForestSetAdaptivityForest()).
112: Collective on dm
114: Input Parameters:
115: + dm - the source DM object
116: - comm - the communicator for the new DM (this communicator is currently ignored, but is present so that DMForestTemplate() can be used within DMCoarsen())
118: Output Parameter:
119: . tdm - the new DM object
121: Level: intermediate
123: .seealso: DMForestSetAdaptivityForest()
124: @*/
125: PetscErrorCode DMForestTemplate(DM dm, MPI_Comm comm, DM *tdm)
126: {
127: DM_Forest *forest = (DM_Forest*) dm->data;
128: DMType type;
129: DM base;
130: DMForestTopology topology;
131: MatType mtype;
132: PetscInt dim, overlap, ref, factor;
133: DMForestAdaptivityStrategy strat;
134: void *ctx;
135: PetscErrorCode (*map)(DM, PetscInt, PetscInt, const PetscReal[], PetscReal[], void*);
136: void *mapCtx;
137: PetscErrorCode ierr;
141: DMCreate(PetscObjectComm((PetscObject)dm),tdm);
142: DMGetType(dm,&type);
143: DMSetType(*tdm,type);
144: DMForestGetBaseDM(dm,&base);
145: DMForestSetBaseDM(*tdm,base);
146: DMForestGetTopology(dm,&topology);
147: DMForestSetTopology(*tdm,topology);
148: DMForestGetAdjacencyDimension(dm,&dim);
149: DMForestSetAdjacencyDimension(*tdm,dim);
150: DMForestGetPartitionOverlap(dm,&overlap);
151: DMForestSetPartitionOverlap(*tdm,overlap);
152: DMForestGetMinimumRefinement(dm,&ref);
153: DMForestSetMinimumRefinement(*tdm,ref);
154: DMForestGetMaximumRefinement(dm,&ref);
155: DMForestSetMaximumRefinement(*tdm,ref);
156: DMForestGetAdaptivityStrategy(dm,&strat);
157: DMForestSetAdaptivityStrategy(*tdm,strat);
158: DMForestGetGradeFactor(dm,&factor);
159: DMForestSetGradeFactor(*tdm,factor);
160: DMForestGetBaseCoordinateMapping(dm,&map,&mapCtx);
161: DMForestSetBaseCoordinateMapping(*tdm,map,mapCtx);
162: if (forest->ftemplate) {
163: (*forest->ftemplate)(dm, *tdm);
164: }
165: DMForestSetAdaptivityForest(*tdm,dm);
166: DMCopyDisc(dm,*tdm);
167: DMGetApplicationContext(dm,&ctx);
168: DMSetApplicationContext(*tdm,&ctx);
169: {
170: PetscBool isper;
171: const PetscReal *maxCell, *L;
172: const DMBoundaryType *bd;
174: DMGetPeriodicity(dm,&isper,&maxCell,&L,&bd);
175: DMSetPeriodicity(*tdm,isper,maxCell,L,bd);
176: }
177: DMCopyBoundary(dm,*tdm);
178: DMGetMatType(dm,&mtype);
179: DMSetMatType(*tdm,mtype);
180: return(0);
181: }
183: static PetscErrorCode DMInitialize_Forest(DM dm);
185: PETSC_EXTERN PetscErrorCode DMClone_Forest(DM dm, DM *newdm)
186: {
187: DM_Forest *forest = (DM_Forest*) dm->data;
188: const char *type;
192: forest->refct++;
193: (*newdm)->data = forest;
194: PetscObjectGetType((PetscObject) dm, &type);
195: PetscObjectChangeTypeName((PetscObject) *newdm, type);
196: DMInitialize_Forest(*newdm);
197: return(0);
198: }
200: static PetscErrorCode DMDestroy_Forest(DM dm)
201: {
202: DM_Forest *forest = (DM_Forest*) dm->data;
206: if (--forest->refct > 0) return(0);
207: if (forest->destroy) {(*forest->destroy)(dm);}
208: PetscSFDestroy(&forest->cellSF);
209: PetscSFDestroy(&forest->preCoarseToFine);
210: PetscSFDestroy(&forest->coarseToPreFine);
211: DMLabelDestroy(&forest->adaptLabel);
212: PetscFree(forest->adaptStrategy);
213: DMDestroy(&forest->base);
214: DMDestroy(&forest->adapt);
215: PetscFree(forest->topology);
216: PetscFree(forest);
217: return(0);
218: }
220: /*@C
221: DMForestSetTopology - Set the topology of a DMForest during the pre-setup phase. The topology is a string (e.g.
222: "cube", "shell") and can be interpreted by subtypes of DMFOREST) to construct the base DM of a forest during
223: DMSetUp().
225: Logically collective on dm
227: Input parameters:
228: + dm - the forest
229: - topology - the topology of the forest
231: Level: intermediate
233: .seealso: DMForestGetTopology(), DMForestSetBaseDM()
234: @*/
235: PetscErrorCode DMForestSetTopology(DM dm, DMForestTopology topology)
236: {
237: DM_Forest *forest = (DM_Forest*) dm->data;
242: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the topology after setup");
243: PetscFree(forest->topology);
244: PetscStrallocpy((const char*)topology,(char**) &forest->topology);
245: return(0);
246: }
248: /*@C
249: DMForestGetTopology - Get a string describing the topology of a DMForest.
251: Not collective
253: Input parameter:
254: . dm - the forest
256: Output parameter:
257: . topology - the topology of the forest (e.g., 'cube', 'shell')
259: Level: intermediate
261: .seealso: DMForestSetTopology()
262: @*/
263: PetscErrorCode DMForestGetTopology(DM dm, DMForestTopology *topology)
264: {
265: DM_Forest *forest = (DM_Forest*) dm->data;
270: *topology = forest->topology;
271: return(0);
272: }
274: /*@
275: DMForestSetBaseDM - During the pre-setup phase, set the DM that defines the base mesh of a DMForest forest. The
276: forest will be hierarchically refined from the base, and all refinements/coarsenings of the forest will share its
277: base. In general, two forest must share a base to be comparable, to do things like construct interpolators.
279: Logically collective on dm
281: Input Parameters:
282: + dm - the forest
283: - base - the base DM of the forest
285: Notes:
286: Currently the base DM must be a DMPLEX
288: Level: intermediate
290: .seealso: DMForestGetBaseDM()
291: @*/
292: PetscErrorCode DMForestSetBaseDM(DM dm, DM base)
293: {
294: DM_Forest *forest = (DM_Forest*) dm->data;
295: PetscInt dim, dimEmbed;
300: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the base after setup");
301: PetscObjectReference((PetscObject)base);
302: DMDestroy(&forest->base);
303: forest->base = base;
304: if (base) {
305: PetscBool isper;
306: const PetscReal *maxCell, *L;
307: const DMBoundaryType *bd;
310: DMGetDimension(base,&dim);
311: DMSetDimension(dm,dim);
312: DMGetCoordinateDim(base,&dimEmbed);
313: DMSetCoordinateDim(dm,dimEmbed);
314: DMGetPeriodicity(base,&isper,&maxCell,&L,&bd);
315: DMSetPeriodicity(dm,isper,maxCell,L,bd);
316: } else {
317: DMSetPeriodicity(dm,PETSC_FALSE,NULL,NULL,NULL);
318: }
319: return(0);
320: }
322: /*@
323: DMForestGetBaseDM - Get the base DM of a DMForest forest. The forest will be hierarchically refined from the base,
324: and all refinements/coarsenings of the forest will share its base. In general, two forest must share a base to be
325: comparable, to do things like construct interpolators.
327: Not collective
329: Input Parameter:
330: . dm - the forest
332: Output Parameter:
333: . base - the base DM of the forest
335: Notes:
336: After DMSetUp(), the base DM will be redundantly distributed across MPI processes
338: Level: intermediate
340: .seealso: DMForestSetBaseDM()
341: @*/
342: PetscErrorCode DMForestGetBaseDM(DM dm, DM *base)
343: {
344: DM_Forest *forest = (DM_Forest*) dm->data;
349: *base = forest->base;
350: return(0);
351: }
353: PetscErrorCode DMForestSetBaseCoordinateMapping(DM dm, PetscErrorCode (*func)(DM,PetscInt,PetscInt,const PetscReal [],PetscReal [],void*),void *ctx)
354: {
355: DM_Forest *forest = (DM_Forest*) dm->data;
359: forest->mapcoordinates = func;
360: forest->mapcoordinatesctx = ctx;
361: return(0);
362: }
364: PetscErrorCode DMForestGetBaseCoordinateMapping(DM dm, PetscErrorCode (**func) (DM,PetscInt,PetscInt,const PetscReal [],PetscReal [],void*),void *ctx)
365: {
366: DM_Forest *forest = (DM_Forest*) dm->data;
370: if (func) *func = forest->mapcoordinates;
371: if (ctx) *((void**) ctx) = forest->mapcoordinatesctx;
372: return(0);
373: }
375: /*@
376: DMForestSetAdaptivityForest - During the pre-setup phase, set the forest from which the current forest will be
377: adapted (e.g., the current forest will be refined/coarsened/repartitioned from it) im DMSetUp(). Usually not needed
378: by users directly: DMForestTemplate() constructs a new forest to be adapted from an old forest and calls this
379: routine.
381: Note that this can be called after setup with adapt = NULL, which will clear all internal data related to the
382: adaptivity forest from dm. This way, repeatedly adapting does not leave stale DM objects in memory.
384: Logically collective on dm
386: Input Parameter:
387: + dm - the new forest, which will be constructed from adapt
388: - adapt - the old forest
390: Level: intermediate
392: .seealso: DMForestGetAdaptivityForest(), DMForestSetAdaptivityPurpose()
393: @*/
394: PetscErrorCode DMForestSetAdaptivityForest(DM dm,DM adapt)
395: {
396: DM_Forest *forest, *adaptForest, *oldAdaptForest;
397: DM oldAdapt;
398: PetscBool isForest;
404: DMIsForest(dm, &isForest);
405: if (!isForest) return(0);
406: if (adapt != NULL && dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the adaptation forest after setup");
407: forest = (DM_Forest*) dm->data;
408: DMForestGetAdaptivityForest(dm,&oldAdapt);
409: adaptForest = (DM_Forest*) (adapt ? adapt->data : NULL);
410: oldAdaptForest = (DM_Forest*) (oldAdapt ? oldAdapt->data : NULL);
411: if (adaptForest != oldAdaptForest) {
412: PetscSFDestroy(&forest->preCoarseToFine);
413: PetscSFDestroy(&forest->coarseToPreFine);
414: if (forest->clearadaptivityforest) {(*forest->clearadaptivityforest)(dm);}
415: }
416: switch (forest->adaptPurpose) {
417: case DM_ADAPT_DETERMINE:
418: PetscObjectReference((PetscObject)adapt);
419: DMDestroy(&(forest->adapt));
420: forest->adapt = adapt;
421: break;
422: case DM_ADAPT_REFINE:
423: DMSetCoarseDM(dm,adapt);
424: break;
425: case DM_ADAPT_COARSEN:
426: case DM_ADAPT_COARSEN_LAST:
427: DMSetFineDM(dm,adapt);
428: break;
429: default:
430: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"invalid adaptivity purpose");
431: }
432: return(0);
433: }
435: /*@
436: DMForestGetAdaptivityForest - Get the forest from which the current forest is adapted.
438: Not collective
440: Input Parameter:
441: . dm - the forest
443: Output Parameter:
444: . adapt - the forest from which dm is/was adapted
446: Level: intermediate
448: .seealso: DMForestSetAdaptivityForest(), DMForestSetAdaptivityPurpose()
449: @*/
450: PetscErrorCode DMForestGetAdaptivityForest(DM dm, DM *adapt)
451: {
452: DM_Forest *forest;
457: forest = (DM_Forest*) dm->data;
458: switch (forest->adaptPurpose) {
459: case DM_ADAPT_DETERMINE:
460: *adapt = forest->adapt;
461: break;
462: case DM_ADAPT_REFINE:
463: DMGetCoarseDM(dm,adapt);
464: break;
465: case DM_ADAPT_COARSEN:
466: case DM_ADAPT_COARSEN_LAST:
467: DMGetFineDM(dm,adapt);
468: break;
469: default:
470: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"invalid adaptivity purpose");
471: }
472: return(0);
473: }
475: /*@
476: DMForestSetAdaptivityPurpose - During the pre-setup phase, set whether the current DM is being adapted from its
477: source (set with DMForestSetAdaptivityForest()) for the purpose of refinement (DM_ADAPT_REFINE), coarsening
478: (DM_ADAPT_COARSEN), or undefined (DM_ADAPT_DETERMINE). This only matters for the purposes of reference counting:
479: during DMDestroy(), cyclic references can be found between DMs only if the cyclic reference is due to a fine/coarse
480: relationship (see DMSetFineDM()/DMSetCoarseDM()). If the purpose is not refinement or coarsening, and the user does
481: not maintain a reference to the post-adaptation forest (i.e., the one created by DMForestTemplate()), then this can
482: cause a memory leak. This method is used by subtypes of DMForest when automatically constructing mesh hierarchies.
484: Logically collective on dm
486: Input Parameters:
487: + dm - the forest
488: - purpose - the adaptivity purpose
490: Level: advanced
492: .seealso: DMForestTemplate(), DMForestSetAdaptivityForest(), DMForestGetAdaptivityForest(), DMAdaptFlag
493: @*/
494: PetscErrorCode DMForestSetAdaptivityPurpose(DM dm, DMAdaptFlag purpose)
495: {
496: DM_Forest *forest;
500: forest = (DM_Forest*) dm->data;
501: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the adaptation forest after setup");
502: if (purpose != forest->adaptPurpose) {
503: DM adapt;
505: DMForestGetAdaptivityForest(dm,&adapt);
506: PetscObjectReference((PetscObject)adapt);
507: DMForestSetAdaptivityForest(dm,NULL);
509: forest->adaptPurpose = purpose;
511: DMForestSetAdaptivityForest(dm,adapt);
512: DMDestroy(&adapt);
513: }
514: return(0);
515: }
517: /*@
518: DMForestGetAdaptivityPurpose - Get whether the current DM is being adapted from its source (set with
519: DMForestSetAdaptivityForest()) for the purpose of refinement (DM_ADAPT_REFINE), coarsening (DM_ADAPT_COARSEN),
520: coarsening only the last level (DM_ADAPT_COARSEN_LAST) or undefined (DM_ADAPT_DETERMINE).
521: This only matters for the purposes of reference counting: during DMDestroy(), cyclic
522: references can be found between DMs only if the cyclic reference is due to a fine/coarse relationship (see
523: DMSetFineDM()/DMSetCoarseDM()). If the purpose is not refinement or coarsening, and the user does not maintain a
524: reference to the post-adaptation forest (i.e., the one created by DMForestTemplate()), then this can cause a memory
525: leak. This method is used by subtypes of DMForest when automatically constructing mesh hierarchies.
527: Not collective
529: Input Parameter:
530: . dm - the forest
532: Output Parameter:
533: . purpose - the adaptivity purpose
535: Level: advanced
537: .seealso: DMForestTemplate(), DMForestSetAdaptivityForest(), DMForestGetAdaptivityForest(), DMAdaptFlag
538: @*/
539: PetscErrorCode DMForestGetAdaptivityPurpose(DM dm, DMAdaptFlag *purpose)
540: {
541: DM_Forest *forest;
544: forest = (DM_Forest*) dm->data;
545: *purpose = forest->adaptPurpose;
546: return(0);
547: }
549: /*@
550: DMForestSetAdjacencyDimension - During the pre-setup phase, set the dimension of interface points that determine
551: cell adjacency (for the purposes of partitioning and overlap).
553: Logically collective on dm
555: Input Parameters:
556: + dm - the forest
557: - adjDim - default 0 (i.e., vertices determine adjacency)
559: Level: intermediate
561: .seealso: DMForestGetAdjacencyDimension(), DMForestSetAdjacencyCodimension(), DMForestSetPartitionOverlap()
562: @*/
563: PetscErrorCode DMForestSetAdjacencyDimension(DM dm, PetscInt adjDim)
564: {
565: PetscInt dim;
566: DM_Forest *forest = (DM_Forest*) dm->data;
571: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the adjacency dimension after setup");
572: if (adjDim < 0) SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_OUTOFRANGE,"adjacency dim cannot be < 0: %d", adjDim);
573: DMGetDimension(dm,&dim);
574: if (adjDim > dim) SETERRQ2(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_OUTOFRANGE,"adjacency dim cannot be > %d: %d", dim, adjDim);
575: forest->adjDim = adjDim;
576: return(0);
577: }
579: /*@
580: DMForestSetAdjacencyCodimension - Like DMForestSetAdjacencyDimension(), but specified as a co-dimension (so that,
581: e.g., adjacency based on facets can be specified by codimension 1 in all cases)
583: Logically collective on dm
585: Input Parameters:
586: + dm - the forest
587: - adjCodim - default isthe dimension of the forest (see DMGetDimension()), since this is the codimension of vertices
589: Level: intermediate
591: .seealso: DMForestGetAdjacencyCodimension(), DMForestSetAdjacencyDimension()
592: @*/
593: PetscErrorCode DMForestSetAdjacencyCodimension(DM dm, PetscInt adjCodim)
594: {
595: PetscInt dim;
600: DMGetDimension(dm,&dim);
601: DMForestSetAdjacencyDimension(dm,dim-adjCodim);
602: return(0);
603: }
605: /*@
606: DMForestGetAdjacencyDimension - Get the dimension of interface points that determine cell adjacency (for the
607: purposes of partitioning and overlap).
609: Not collective
611: Input Parameter:
612: . dm - the forest
614: Output Parameter:
615: . adjDim - default 0 (i.e., vertices determine adjacency)
617: Level: intermediate
619: .seealso: DMForestSetAdjacencyDimension(), DMForestGetAdjacencyCodimension(), DMForestSetPartitionOverlap()
620: @*/
621: PetscErrorCode DMForestGetAdjacencyDimension(DM dm, PetscInt *adjDim)
622: {
623: DM_Forest *forest = (DM_Forest*) dm->data;
628: *adjDim = forest->adjDim;
629: return(0);
630: }
632: /*@
633: DMForestGetAdjacencyCodimension - Like DMForestGetAdjacencyDimension(), but specified as a co-dimension (so that,
634: e.g., adjacency based on facets can be specified by codimension 1 in all cases)
636: Not collective
638: Input Parameter:
639: . dm - the forest
641: Output Parameter:
642: . adjCodim - default isthe dimension of the forest (see DMGetDimension()), since this is the codimension of vertices
644: Level: intermediate
646: .seealso: DMForestSetAdjacencyCodimension(), DMForestGetAdjacencyDimension()
647: @*/
648: PetscErrorCode DMForestGetAdjacencyCodimension(DM dm, PetscInt *adjCodim)
649: {
650: DM_Forest *forest = (DM_Forest*) dm->data;
651: PetscInt dim;
657: DMGetDimension(dm,&dim);
658: *adjCodim = dim - forest->adjDim;
659: return(0);
660: }
662: /*@
663: DMForestSetPartitionOverlap - During the pre-setup phase, set the amount of cell-overlap present in parallel
664: partitions of a forest, with values > 0 indicating subdomains that are expanded by that many iterations of adding
665: adjacent cells
667: Logically collective on dm
669: Input Parameters:
670: + dm - the forest
671: - overlap - default 0
673: Level: intermediate
675: .seealso: DMForestGetPartitionOverlap(), DMForestSetAdjacencyDimension(), DMForestSetAdjacencyCodimension()
676: @*/
677: PetscErrorCode DMForestSetPartitionOverlap(DM dm, PetscInt overlap)
678: {
679: DM_Forest *forest = (DM_Forest*) dm->data;
683: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the overlap after setup");
684: if (overlap < 0) SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_OUTOFRANGE,"overlap cannot be < 0: %d", overlap);
685: forest->overlap = overlap;
686: return(0);
687: }
689: /*@
690: DMForestGetPartitionOverlap - Get the amount of cell-overlap present in parallel partitions of a forest, with values
691: > 0 indicating subdomains that are expanded by that many iterations of adding adjacent cells
693: Not collective
695: Input Parameter:
696: . dm - the forest
698: Output Parameter:
699: . overlap - default 0
701: Level: intermediate
703: .seealso: DMForestGetPartitionOverlap(), DMForestSetAdjacencyDimension(), DMForestSetAdjacencyCodimension()
704: @*/
705: PetscErrorCode DMForestGetPartitionOverlap(DM dm, PetscInt *overlap)
706: {
707: DM_Forest *forest = (DM_Forest*) dm->data;
712: *overlap = forest->overlap;
713: return(0);
714: }
716: /*@
717: DMForestSetMinimumRefinement - During the pre-setup phase, set the minimum level of refinement (relative to the base
718: DM, see DMForestGetBaseDM()) allowed in the forest. If the forest is being created by coarsening a previous forest
719: (see DMForestGetAdaptivityForest()) this limits the amount of coarsening.
721: Logically collective on dm
723: Input Parameters:
724: + dm - the forest
725: - minRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
727: Level: intermediate
729: .seealso: DMForestGetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestSetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityForest()
730: @*/
731: PetscErrorCode DMForestSetMinimumRefinement(DM dm, PetscInt minRefinement)
732: {
733: DM_Forest *forest = (DM_Forest*) dm->data;
737: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the minimum refinement after setup");
738: forest->minRefinement = minRefinement;
739: return(0);
740: }
742: /*@
743: DMForestGetMinimumRefinement - Get the minimum level of refinement (relative to the base DM, see
744: DMForestGetBaseDM()) allowed in the forest. If the forest is being created by coarsening a previous forest (see
745: DMForestGetAdaptivityForest()), this limits the amount of coarsening.
747: Not collective
749: Input Parameter:
750: . dm - the forest
752: Output Parameter:
753: . minRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
755: Level: intermediate
757: .seealso: DMForestSetMinimumRefinement(), DMForestGetMaximumRefinement(), DMForestGetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityForest()
758: @*/
759: PetscErrorCode DMForestGetMinimumRefinement(DM dm, PetscInt *minRefinement)
760: {
761: DM_Forest *forest = (DM_Forest*) dm->data;
766: *minRefinement = forest->minRefinement;
767: return(0);
768: }
770: /*@
771: DMForestSetInitialRefinement - During the pre-setup phase, set the initial level of refinement (relative to the base
772: DM, see DMForestGetBaseDM()) allowed in the forest.
774: Logically collective on dm
776: Input Parameters:
777: + dm - the forest
778: - initefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
780: Level: intermediate
782: .seealso: DMForestSetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestGetBaseDM()
783: @*/
784: PetscErrorCode DMForestSetInitialRefinement(DM dm, PetscInt initRefinement)
785: {
786: DM_Forest *forest = (DM_Forest*) dm->data;
790: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the initial refinement after setup");
791: forest->initRefinement = initRefinement;
792: return(0);
793: }
795: /*@
796: DMForestGetInitialRefinement - Get the initial level of refinement (relative to the base DM, see
797: DMForestGetBaseDM()) allowed in the forest.
799: Not collective
801: Input Parameter:
802: . dm - the forest
804: Output Paramater:
805: . initRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
807: Level: intermediate
809: .seealso: DMForestSetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestGetBaseDM()
810: @*/
811: PetscErrorCode DMForestGetInitialRefinement(DM dm, PetscInt *initRefinement)
812: {
813: DM_Forest *forest = (DM_Forest*) dm->data;
818: *initRefinement = forest->initRefinement;
819: return(0);
820: }
822: /*@
823: DMForestSetMaximumRefinement - During the pre-setup phase, set the maximum level of refinement (relative to the base
824: DM, see DMForestGetBaseDM()) allowed in the forest. If the forest is being created by refining a previous forest
825: (see DMForestGetAdaptivityForest()), this limits the amount of refinement.
827: Logically collective on dm
829: Input Parameters:
830: + dm - the forest
831: - maxRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
833: Level: intermediate
835: .seealso: DMForestGetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestSetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityDM()
836: @*/
837: PetscErrorCode DMForestSetMaximumRefinement(DM dm, PetscInt maxRefinement)
838: {
839: DM_Forest *forest = (DM_Forest*) dm->data;
843: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the maximum refinement after setup");
844: forest->maxRefinement = maxRefinement;
845: return(0);
846: }
848: /*@
849: DMForestGetMaximumRefinement - Get the maximum level of refinement (relative to the base DM, see
850: DMForestGetBaseDM()) allowed in the forest. If the forest is being created by refining a previous forest (see
851: DMForestGetAdaptivityForest()), this limits the amount of refinement.
853: Not collective
855: Input Parameter:
856: . dm - the forest
858: Output Parameter:
859: . maxRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
861: Level: intermediate
863: .seealso: DMForestSetMaximumRefinement(), DMForestGetMinimumRefinement(), DMForestGetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityForest()
864: @*/
865: PetscErrorCode DMForestGetMaximumRefinement(DM dm, PetscInt *maxRefinement)
866: {
867: DM_Forest *forest = (DM_Forest*) dm->data;
872: *maxRefinement = forest->maxRefinement;
873: return(0);
874: }
876: /*@C
877: DMForestSetAdaptivityStrategy - During the pre-setup phase, set the strategy for combining adaptivity labels from multiple processes.
878: Subtypes of DMForest may define their own strategies. Two default strategies are DMFORESTADAPTALL, which indicates that all processes must agree
879: for a refinement/coarsening flag to be valid, and DMFORESTADAPTANY, which indicates that only one process needs to
880: specify refinement/coarsening.
882: Logically collective on dm
884: Input Parameters:
885: + dm - the forest
886: - adaptStrategy - default DMFORESTADAPTALL
888: Level: advanced
890: .seealso: DMForestGetAdaptivityStrategy()
891: @*/
892: PetscErrorCode DMForestSetAdaptivityStrategy(DM dm, DMForestAdaptivityStrategy adaptStrategy)
893: {
894: DM_Forest *forest = (DM_Forest*) dm->data;
899: PetscFree(forest->adaptStrategy);
900: PetscStrallocpy((const char*) adaptStrategy,(char**)&forest->adaptStrategy);
901: return(0);
902: }
904: /*@C
905: DMForestSetAdaptivityStrategy - Get the strategy for combining adaptivity labels from multiple processes. Subtypes
906: of DMForest may define their own strategies. Two default strategies are DMFORESTADAPTALL, which indicates that all
907: processes must agree for a refinement/coarsening flag to be valid, and DMFORESTADAPTANY, which indicates that only
908: one process needs to specify refinement/coarsening.
910: Not collective
912: Input Parameter:
913: . dm - the forest
915: Output Parameter:
916: . adaptStrategy - the adaptivity strategy (default DMFORESTADAPTALL)
918: Level: advanced
920: .seealso: DMForestSetAdaptivityStrategy()
921: @*/
922: PetscErrorCode DMForestGetAdaptivityStrategy(DM dm, DMForestAdaptivityStrategy *adaptStrategy)
923: {
924: DM_Forest *forest = (DM_Forest*) dm->data;
929: *adaptStrategy = forest->adaptStrategy;
930: return(0);
931: }
933: /*@
934: DMForestGetAdaptivitySuccess - Return whether the requested adaptation (refinement, coarsening, repartitioning,
935: etc.) was successful. PETSC_FALSE indicates that the post-adaptation forest is the same as the pre-adpatation
936: forest. A requested adaptation may have been unsuccessful if, for example, the requested refinement would have
937: exceeded the maximum refinement level.
939: Collective on dm
941: Input Parameter:
943: . dm - the post-adaptation forest
945: Output Parameter:
947: . success - PETSC_TRUE if the post-adaptation forest is different from the pre-adaptation forest.
949: Level: intermediate
951: .see
952: @*/
953: PetscErrorCode DMForestGetAdaptivitySuccess(DM dm, PetscBool *success)
954: {
955: DM_Forest *forest;
960: if (!dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"DMSetUp() has not been called yet.");
961: forest = (DM_Forest *) dm->data;
962: (forest->getadaptivitysuccess)(dm,success);
963: return(0);
964: }
966: /*@
967: DMForestSetComputeAdaptivitySF - During the pre-setup phase, set whether transfer PetscSFs should be computed
968: 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().
970: Logically collective on dm
972: Input Parameters:
973: + dm - the post-adaptation forest
974: - computeSF - default PETSC_TRUE
976: Level: advanced
978: .seealso: DMForestGetComputeAdaptivitySF(), DMForestGetAdaptivitySF()
979: @*/
980: PetscErrorCode DMForestSetComputeAdaptivitySF(DM dm, PetscBool computeSF)
981: {
982: DM_Forest *forest;
986: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot compute adaptivity PetscSFs after setup is called");
987: forest = (DM_Forest*) dm->data;
988: forest->computeAdaptSF = computeSF;
989: return(0);
990: }
992: PetscErrorCode DMForestTransferVec(DM dmIn, Vec vecIn, DM dmOut, Vec vecOut, PetscBool useBCs, PetscReal time)
993: {
994: DM_Forest *forest;
1002: forest = (DM_Forest *) dmIn->data;
1003: if (!forest->transfervec) SETERRQ(PetscObjectComm((PetscObject)dmIn),PETSC_ERR_SUP,"DMForestTransferVec() not implemented");
1004: (forest->transfervec)(dmIn,vecIn,dmOut,vecOut,useBCs,time);
1005: return(0);
1006: }
1008: PetscErrorCode DMForestTransferVecFromBase(DM dm, Vec vecIn, Vec vecOut)
1009: {
1010: DM_Forest *forest;
1017: forest = (DM_Forest *) dm->data;
1018: if (!forest->transfervecfrombase) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"DMForestTransferVecFromBase() not implemented");
1019: (forest->transfervecfrombase)(dm,vecIn,vecOut);
1020: return(0);
1021: }
1023: /*@
1024: DMForestGetComputeAdaptivitySF - Get whether transfer PetscSFs should be computed relating the cells of the
1025: pre-adaptation forest to the post-adaptiation forest. After DMSetUp() is called, these transfer PetscSFs can be
1026: accessed with DMForestGetAdaptivitySF().
1028: Not collective
1030: Input Parameter:
1031: . dm - the post-adaptation forest
1033: Output Parameter:
1034: . computeSF - default PETSC_TRUE
1036: Level: advanced
1038: .seealso: DMForestSetComputeAdaptivitySF(), DMForestGetAdaptivitySF()
1039: @*/
1040: PetscErrorCode DMForestGetComputeAdaptivitySF(DM dm, PetscBool *computeSF)
1041: {
1042: DM_Forest *forest;
1046: forest = (DM_Forest*) dm->data;
1047: *computeSF = forest->computeAdaptSF;
1048: return(0);
1049: }
1051: /*@
1052: DMForestGetAdaptivitySF - Get PetscSFs that relate the pre-adaptation forest to the post-adaptation forest.
1053: Adaptation can be any combination of refinement, coarsening, repartition, and change of overlap, so there may be
1054: some cells of the pre-adaptation that are parents of post-adaptation cells, and vice versa. Therefore there are two
1055: PetscSFs: one that relates pre-adaptation coarse cells to post-adaptation fine cells, and one that relates
1056: pre-adaptation fine cells to post-adaptation coarse cells.
1058: Not collective
1060: Input Parameter:
1061: dm - the post-adaptation forest
1063: Output Parameter:
1064: preCoarseToFine - pre-adaptation coarse cells to post-adaptation fine cells: BCast goes from pre- to post-
1065: coarseToPreFine - post-adaptation coarse cells to pre-adaptation fine cells: BCast goes from post- to pre-
1067: Level: advanced
1069: .seealso: DMForestGetComputeAdaptivitySF(), DMForestSetComputeAdaptivitySF()
1070: @*/
1071: PetscErrorCode DMForestGetAdaptivitySF(DM dm, PetscSF *preCoarseToFine, PetscSF *coarseToPreFine)
1072: {
1073: DM_Forest *forest;
1078: DMSetUp(dm);
1079: forest = (DM_Forest*) dm->data;
1080: if (preCoarseToFine) *preCoarseToFine = forest->preCoarseToFine;
1081: if (coarseToPreFine) *coarseToPreFine = forest->coarseToPreFine;
1082: return(0);
1083: }
1085: /*@
1086: DMForestSetGradeFactor - During the pre-setup phase, set the desired amount of grading in the mesh, e.g. give 2 to
1087: indicate that the diameter of neighboring cells should differ by at most a factor of 2. Subtypes of DMForest may
1088: only support one particular choice of grading factor.
1090: Logically collective on dm
1092: Input Parameters:
1093: + dm - the forest
1094: - grade - the grading factor
1096: Level: advanced
1098: .seealso: DMForestGetGradeFactor()
1099: @*/
1100: PetscErrorCode DMForestSetGradeFactor(DM dm, PetscInt grade)
1101: {
1102: DM_Forest *forest = (DM_Forest*) dm->data;
1106: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the grade factor after setup");
1107: forest->gradeFactor = grade;
1108: return(0);
1109: }
1111: /*@
1112: DMForestGetGradeFactor - Get the desired amount of grading in the mesh, e.g. give 2 to indicate that the diameter of
1113: neighboring cells should differ by at most a factor of 2. Subtypes of DMForest may only support one particular
1114: choice of grading factor.
1116: Not collective
1118: Input Parameter:
1119: . dm - the forest
1121: Output Parameter:
1122: . grade - the grading factor
1124: Level: advanced
1126: .seealso: DMForestSetGradeFactor()
1127: @*/
1128: PetscErrorCode DMForestGetGradeFactor(DM dm, PetscInt *grade)
1129: {
1130: DM_Forest *forest = (DM_Forest*) dm->data;
1135: *grade = forest->gradeFactor;
1136: return(0);
1137: }
1139: /*@
1140: DMForestSetCellWeightFactor - During the pre-setup phase, set the factor by which the level of refinement changes
1141: the cell weight (see DMForestSetCellWeights()) when calculating partitions. The final weight of a cell will be
1142: (cellWeight) * (weightFactor^refinementLevel). A factor of 1 indicates that the weight of a cell does not depend on
1143: its level; a factor of 2, for example, might be appropriate for sub-cycling time-stepping methods, when the
1144: computation associated with a cell is multiplied by a factor of 2 for each additional level of refinement.
1146: Logically collective on dm
1148: Input Parameters:
1149: + dm - the forest
1150: - weightsFactors - default 1.
1152: Level: advanced
1154: .seealso: DMForestGetCellWeightFactor(), DMForestSetCellWeights()
1155: @*/
1156: PetscErrorCode DMForestSetCellWeightFactor(DM dm, PetscReal weightsFactor)
1157: {
1158: DM_Forest *forest = (DM_Forest*) dm->data;
1162: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the weights factor after setup");
1163: forest->weightsFactor = weightsFactor;
1164: return(0);
1165: }
1167: /*@
1168: DMForestGetCellWeightFactor - Get the factor by which the level of refinement changes the cell weight (see
1169: DMForestSetCellWeights()) when calculating partitions. The final weight of a cell will be (cellWeight) *
1170: (weightFactor^refinementLevel). A factor of 1 indicates that the weight of a cell does not depend on its level; a
1171: factor of 2, for example, might be appropriate for sub-cycling time-stepping methods, when the computation
1172: associated with a cell is multiplied by a factor of 2 for each additional level of refinement.
1174: Not collective
1176: Input Parameter:
1177: . dm - the forest
1179: Output Parameter:
1180: . weightsFactors - default 1.
1182: Level: advanced
1184: .seealso: DMForestSetCellWeightFactor(), DMForestSetCellWeights()
1185: @*/
1186: PetscErrorCode DMForestGetCellWeightFactor(DM dm, PetscReal *weightsFactor)
1187: {
1188: DM_Forest *forest = (DM_Forest*) dm->data;
1193: *weightsFactor = forest->weightsFactor;
1194: return(0);
1195: }
1197: /*@
1198: DMForestGetCellChart - After the setup phase, get the local half-open interval of the chart of cells on this process
1200: Not collective
1202: Input Parameter:
1203: . dm - the forest
1205: Output Parameters:
1206: + cStart - the first cell on this process
1207: - cEnd - one after the final cell on this process
1209: Level: intermediate
1211: .seealso: DMForestGetCellSF()
1212: @*/
1213: PetscErrorCode DMForestGetCellChart(DM dm, PetscInt *cStart, PetscInt *cEnd)
1214: {
1215: DM_Forest *forest = (DM_Forest*) dm->data;
1222: if (((forest->cStart == PETSC_DETERMINE) || (forest->cEnd == PETSC_DETERMINE)) && forest->createcellchart) {
1223: forest->createcellchart(dm,&forest->cStart,&forest->cEnd);
1224: }
1225: *cStart = forest->cStart;
1226: *cEnd = forest->cEnd;
1227: return(0);
1228: }
1230: /*@
1231: DMForestGetCellSF - After the setup phase, get the PetscSF for overlapping cells between processes
1233: Not collective
1235: Input Parameter:
1236: . dm - the forest
1238: Output Parameter:
1239: . cellSF - the PetscSF
1241: Level: intermediate
1243: .seealso: DMForestGetCellChart()
1244: @*/
1245: PetscErrorCode DMForestGetCellSF(DM dm, PetscSF *cellSF)
1246: {
1247: DM_Forest *forest = (DM_Forest*) dm->data;
1253: if ((!forest->cellSF) && forest->createcellsf) {
1254: forest->createcellsf(dm,&forest->cellSF);
1255: }
1256: *cellSF = forest->cellSF;
1257: return(0);
1258: }
1260: /*@C
1261: DMForestSetAdaptivityLabel - During the pre-setup phase, set the label of the pre-adaptation forest (see
1262: DMForestGetAdaptivityForest()) that holds the adaptation flags (refinement, coarsening, or some combination). The
1263: interpretation of the label values is up to the subtype of DMForest, but DM_ADAPT_DETERMINE, DM_ADAPT_KEEP,
1264: DM_ADAPT_REFINE, and DM_ADAPT_COARSEN have been reserved as choices that should be accepted by all subtypes.
1266: Logically collective on dm
1268: Input Parameters:
1269: - dm - the forest
1270: + adaptLabel - the label in the pre-adaptation forest
1272: Level: intermediate
1274: .seealso DMForestGetAdaptivityLabel()
1275: @*/
1276: PetscErrorCode DMForestSetAdaptivityLabel(DM dm, DMLabel adaptLabel)
1277: {
1278: DM_Forest *forest = (DM_Forest*) dm->data;
1284: PetscObjectReference((PetscObject)adaptLabel);
1285: DMLabelDestroy(&forest->adaptLabel);
1286: forest->adaptLabel = adaptLabel;
1287: return(0);
1288: }
1290: /*@C
1291: DMForestGetAdaptivityLabel - Get the label of the pre-adaptation forest (see DMForestGetAdaptivityForest()) that
1292: holds the adaptation flags (refinement, coarsening, or some combination). The interpretation of the label values is
1293: up to the subtype of DMForest, but DM_ADAPT_DETERMINE, DM_ADAPT_KEEP, DM_ADAPT_REFINE, and DM_ADAPT_COARSEN have
1294: been reserved as choices that should be accepted by all subtypes.
1296: Not collective
1298: Input Parameter:
1299: . dm - the forest
1301: Output Parameter:
1302: . adaptLabel - the name of the label in the pre-adaptation forest
1304: Level: intermediate
1306: .seealso DMForestSetAdaptivityLabel()
1307: @*/
1308: PetscErrorCode DMForestGetAdaptivityLabel(DM dm, DMLabel *adaptLabel)
1309: {
1310: DM_Forest *forest = (DM_Forest*) dm->data;
1314: *adaptLabel = forest->adaptLabel;
1315: return(0);
1316: }
1318: /*@
1319: DMForestSetCellWeights - Set 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: Logically collective on dm
1326: Input Parameters:
1327: + dm - the forest
1328: . weights - the array of weights for all cells, or NULL to indicate each cell has weight 1.
1329: - copyMode - how weights should reference weights
1331: Level: advanced
1333: .seealso: DMForestGetCellWeights(), DMForestSetWeightCapacity()
1334: @*/
1335: PetscErrorCode DMForestSetCellWeights(DM dm, PetscReal weights[], PetscCopyMode copyMode)
1336: {
1337: DM_Forest *forest = (DM_Forest*) dm->data;
1338: PetscInt cStart, cEnd;
1343: DMForestGetCellChart(dm,&cStart,&cEnd);
1344: if (cEnd < cStart) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"cell chart [%d,%d) is not valid",cStart,cEnd);
1345: if (copyMode == PETSC_COPY_VALUES) {
1346: if (forest->cellWeightsCopyMode != PETSC_OWN_POINTER || forest->cellWeights == weights) {
1347: PetscMalloc1(cEnd-cStart,&forest->cellWeights);
1348: }
1349: PetscArraycpy(forest->cellWeights,weights,cEnd-cStart);
1350: forest->cellWeightsCopyMode = PETSC_OWN_POINTER;
1351: return(0);
1352: }
1353: if (forest->cellWeightsCopyMode == PETSC_OWN_POINTER) {
1354: PetscFree(forest->cellWeights);
1355: }
1356: forest->cellWeights = weights;
1357: forest->cellWeightsCopyMode = copyMode;
1358: return(0);
1359: }
1361: /*@
1362: DMForestGetCellWeights - Get the weights assigned to each of the cells (see DMForestGetCellChart()) of the current
1363: process: weights are used to determine parallel partitioning. Partitions will be created so that each process's
1364: ratio of weight to capacity (see DMForestSetWeightCapacity()) is roughly equal. If NULL, each cell receives a weight
1365: of 1.
1367: Not collective
1369: Input Parameter:
1370: . dm - the forest
1372: Output Parameter:
1373: . weights - the array of weights for all cells, or NULL to indicate each cell has weight 1.
1375: Level: advanced
1377: .seealso: DMForestSetCellWeights(), DMForestSetWeightCapacity()
1378: @*/
1379: PetscErrorCode DMForestGetCellWeights(DM dm, PetscReal **weights)
1380: {
1381: DM_Forest *forest = (DM_Forest*) dm->data;
1386: *weights = forest->cellWeights;
1387: return(0);
1388: }
1390: /*@
1391: DMForestSetWeightCapacity - During the pre-setup phase, set the capacity of the current process when repartitioning
1392: a pre-adaptation forest (see DMForestGetAdaptivityForest()). After partitioning, the ratio of the weight of each
1393: process's cells to the process's capacity will be roughly equal for all processes. A capacity of 0 indicates that
1394: the current process should not have any cells after repartitioning.
1396: Logically Collective on dm
1398: Input parameters:
1399: + dm - the forest
1400: - capacity - this process's capacity
1402: Level: advanced
1404: .seealso DMForestGetWeightCapacity(), DMForestSetCellWeights(), DMForestSetCellWeightFactor()
1405: @*/
1406: PetscErrorCode DMForestSetWeightCapacity(DM dm, PetscReal capacity)
1407: {
1408: DM_Forest *forest = (DM_Forest*) dm->data;
1412: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the weight capacity after setup");
1413: if (capacity < 0.) SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_OUTOFRANGE,"Cannot have negative weight capacity; %f",capacity);
1414: forest->weightCapacity = capacity;
1415: return(0);
1416: }
1418: /*@
1419: DMForestGetWeightCapacity - Set the capacity of the current process when repartitioning a pre-adaptation forest (see
1420: DMForestGetAdaptivityForest()). After partitioning, the ratio of the weight of each process's cells to the
1421: process's capacity will be roughly equal for all processes. A capacity of 0 indicates that the current process
1422: should not have any cells after repartitioning.
1424: Not collective
1426: Input parameter:
1427: . dm - the forest
1429: Output parameter:
1430: . capacity - this process's capacity
1432: Level: advanced
1434: .seealso DMForestSetWeightCapacity(), DMForestSetCellWeights(), DMForestSetCellWeightFactor()
1435: @*/
1436: PetscErrorCode DMForestGetWeightCapacity(DM dm, PetscReal *capacity)
1437: {
1438: DM_Forest *forest = (DM_Forest*) dm->data;
1443: *capacity = forest->weightCapacity;
1444: return(0);
1445: }
1447: PETSC_EXTERN PetscErrorCode DMSetFromOptions_Forest(PetscOptionItems *PetscOptionsObject,DM dm)
1448: {
1449: PetscBool flg, flg1, flg2, flg3, flg4;
1450: DMForestTopology oldTopo;
1451: char stringBuffer[256];
1452: PetscViewer viewer;
1453: PetscViewerFormat format;
1454: PetscInt adjDim, adjCodim, overlap, minRefinement, initRefinement, maxRefinement, grade;
1455: PetscReal weightsFactor;
1456: DMForestAdaptivityStrategy adaptStrategy;
1457: PetscErrorCode ierr;
1461: DMForestGetTopology(dm, &oldTopo);
1462: PetscOptionsHead(PetscOptionsObject,"DMForest Options");
1463: PetscOptionsString("-dm_forest_topology","the topology of the forest's base mesh","DMForestSetTopology",oldTopo,stringBuffer,sizeof(stringBuffer),&flg1);
1464: PetscOptionsViewer("-dm_forest_base_dm","load the base DM from a viewer specification","DMForestSetBaseDM",&viewer,&format,&flg2);
1465: PetscOptionsViewer("-dm_forest_coarse_forest","load the coarse forest from a viewer specification","DMForestSetCoarseForest",&viewer,&format,&flg3);
1466: PetscOptionsViewer("-dm_forest_fine_forest","load the fine forest from a viewer specification","DMForestSetFineForest",&viewer,&format,&flg4);
1467: if ((PetscInt) flg1 + (PetscInt) flg2 + (PetscInt) flg3 + (PetscInt) flg4 > 1) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_INCOMP,"Specify only one of -dm_forest_{topology,base_dm,coarse_forest,fine_forest}");
1468: if (flg1) {
1469: DMForestSetTopology(dm,(DMForestTopology)stringBuffer);
1470: DMForestSetBaseDM(dm,NULL);
1471: DMForestSetAdaptivityForest(dm,NULL);
1472: }
1473: if (flg2) {
1474: DM base;
1476: DMCreate(PetscObjectComm((PetscObject)dm),&base);
1477: PetscViewerPushFormat(viewer,format);
1478: DMLoad(base,viewer);
1479: PetscViewerDestroy(&viewer);
1480: DMForestSetBaseDM(dm,base);
1481: DMDestroy(&base);
1482: DMForestSetTopology(dm,NULL);
1483: DMForestSetAdaptivityForest(dm,NULL);
1484: }
1485: if (flg3) {
1486: DM coarse;
1488: DMCreate(PetscObjectComm((PetscObject)dm),&coarse);
1489: PetscViewerPushFormat(viewer,format);
1490: DMLoad(coarse,viewer);
1491: PetscViewerDestroy(&viewer);
1492: DMForestSetAdaptivityForest(dm,coarse);
1493: DMDestroy(&coarse);
1494: DMForestSetTopology(dm,NULL);
1495: DMForestSetBaseDM(dm,NULL);
1496: }
1497: if (flg4) {
1498: DM fine;
1500: DMCreate(PetscObjectComm((PetscObject)dm),&fine);
1501: PetscViewerPushFormat(viewer,format);
1502: DMLoad(fine,viewer);
1503: PetscViewerDestroy(&viewer);
1504: DMForestSetAdaptivityForest(dm,fine);
1505: DMDestroy(&fine);
1506: DMForestSetTopology(dm,NULL);
1507: DMForestSetBaseDM(dm,NULL);
1508: }
1509: DMForestGetAdjacencyDimension(dm,&adjDim);
1510: PetscOptionsBoundedInt("-dm_forest_adjacency_dimension","set the dimension of points that define adjacency in the forest","DMForestSetAdjacencyDimension",adjDim,&adjDim,&flg,0);
1511: if (flg) {
1512: DMForestSetAdjacencyDimension(dm,adjDim);
1513: } else {
1514: DMForestGetAdjacencyCodimension(dm,&adjCodim);
1515: PetscOptionsBoundedInt("-dm_forest_adjacency_codimension","set the codimension of points that define adjacency in the forest","DMForestSetAdjacencyCodimension",adjCodim,&adjCodim,&flg,1);
1516: if (flg) {
1517: DMForestSetAdjacencyCodimension(dm,adjCodim);
1518: }
1519: }
1520: DMForestGetPartitionOverlap(dm,&overlap);
1521: PetscOptionsBoundedInt("-dm_forest_partition_overlap","set the degree of partition overlap","DMForestSetPartitionOverlap",overlap,&overlap,&flg,0);
1522: if (flg) {
1523: DMForestSetPartitionOverlap(dm,overlap);
1524: }
1525: #if 0
1526: 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);
1527: if (flg) {
1528: DMForestSetMinimumRefinement(dm,minRefinement);
1529: DMForestSetInitialRefinement(dm,minRefinement);
1530: }
1531: PetscOptionsBoundedInt("-dm_refine_hierarchy","equivalent to -dm_forest_set_minimum_refinement 0 and -dm_forest_set_initial_refinement",NULL,initRefinement,&initRefinement,&flg,0);
1532: if (flg) {
1533: DMForestSetMinimumRefinement(dm,0);
1534: DMForestSetInitialRefinement(dm,initRefinement);
1535: }
1536: #endif
1537: DMForestGetMinimumRefinement(dm,&minRefinement);
1538: PetscOptionsBoundedInt("-dm_forest_minimum_refinement","set the minimum level of refinement in the forest","DMForestSetMinimumRefinement",minRefinement,&minRefinement,&flg,0);
1539: if (flg) {
1540: DMForestSetMinimumRefinement(dm,minRefinement);
1541: }
1542: DMForestGetInitialRefinement(dm,&initRefinement);
1543: PetscOptionsBoundedInt("-dm_forest_initial_refinement","set the initial level of refinement in the forest","DMForestSetInitialRefinement",initRefinement,&initRefinement,&flg,0);
1544: if (flg) {
1545: DMForestSetInitialRefinement(dm,initRefinement);
1546: }
1547: DMForestGetMaximumRefinement(dm,&maxRefinement);
1548: PetscOptionsBoundedInt("-dm_forest_maximum_refinement","set the maximum level of refinement in the forest","DMForestSetMaximumRefinement",maxRefinement,&maxRefinement,&flg,0);
1549: if (flg) {
1550: DMForestSetMaximumRefinement(dm,maxRefinement);
1551: }
1552: DMForestGetAdaptivityStrategy(dm,&adaptStrategy);
1553: PetscOptionsString("-dm_forest_adaptivity_strategy","the forest's adaptivity-flag resolution strategy","DMForestSetAdaptivityStrategy",adaptStrategy,stringBuffer,sizeof(stringBuffer),&flg);
1554: if (flg) {
1555: DMForestSetAdaptivityStrategy(dm,(DMForestAdaptivityStrategy)stringBuffer);
1556: }
1557: DMForestGetGradeFactor(dm,&grade);
1558: PetscOptionsBoundedInt("-dm_forest_grade_factor","grade factor between neighboring cells","DMForestSetGradeFactor",grade,&grade,&flg,0);
1559: if (flg) {
1560: DMForestSetGradeFactor(dm,grade);
1561: }
1562: DMForestGetCellWeightFactor(dm,&weightsFactor);
1563: PetscOptionsReal("-dm_forest_cell_weight_factor","multiplying weight factor for cell refinement","DMForestSetCellWeightFactor",weightsFactor,&weightsFactor,&flg);
1564: if (flg) {
1565: DMForestSetCellWeightFactor(dm,weightsFactor);
1566: }
1567: PetscOptionsTail();
1568: return(0);
1569: }
1571: PetscErrorCode DMCreateSubDM_Forest(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1572: {
1576: if (subdm) {DMClone(dm, subdm);}
1577: DMCreateSectionSubDM(dm, numFields, fields, is, subdm);
1578: return(0);
1579: }
1581: PetscErrorCode DMRefine_Forest(DM dm, MPI_Comm comm, DM *dmRefined)
1582: {
1583: DMLabel refine;
1584: DM fineDM;
1588: DMGetFineDM(dm,&fineDM);
1589: if (fineDM) {
1590: PetscObjectReference((PetscObject)fineDM);
1591: *dmRefined = fineDM;
1592: return(0);
1593: }
1594: DMForestTemplate(dm,comm,dmRefined);
1595: DMGetLabel(dm,"refine",&refine);
1596: if (!refine) {
1597: DMLabelCreate(PETSC_COMM_SELF, "refine",&refine);
1598: DMLabelSetDefaultValue(refine,DM_ADAPT_REFINE);
1599: } else {
1600: PetscObjectReference((PetscObject) refine);
1601: }
1602: DMForestSetAdaptivityLabel(*dmRefined,refine);
1603: DMLabelDestroy(&refine);
1604: return(0);
1605: }
1607: PetscErrorCode DMCoarsen_Forest(DM dm, MPI_Comm comm, DM *dmCoarsened)
1608: {
1609: DMLabel coarsen;
1610: DM coarseDM;
1614: {
1615: PetscMPIInt mpiComparison;
1616: MPI_Comm dmcomm = PetscObjectComm((PetscObject)dm);
1618: MPI_Comm_compare(comm, dmcomm, &mpiComparison);
1619: if (mpiComparison != MPI_IDENT && mpiComparison != MPI_CONGRUENT) SETERRQ(dmcomm,PETSC_ERR_SUP,"No support for different communicators yet");
1620: }
1621: DMGetCoarseDM(dm,&coarseDM);
1622: if (coarseDM) {
1623: PetscObjectReference((PetscObject)coarseDM);
1624: *dmCoarsened = coarseDM;
1625: return(0);
1626: }
1627: DMForestTemplate(dm,comm,dmCoarsened);
1628: DMForestSetAdaptivityPurpose(*dmCoarsened,DM_ADAPT_COARSEN);
1629: DMGetLabel(dm,"coarsen",&coarsen);
1630: if (!coarsen) {
1631: DMLabelCreate(PETSC_COMM_SELF, "coarsen",&coarsen);
1632: DMLabelSetDefaultValue(coarsen,DM_ADAPT_COARSEN);
1633: } else {
1634: PetscObjectReference((PetscObject) coarsen);
1635: }
1636: DMForestSetAdaptivityLabel(*dmCoarsened,coarsen);
1637: DMLabelDestroy(&coarsen);
1638: return(0);
1639: }
1641: static PetscErrorCode DMAdaptLabel_Forest(DM dm, DMLabel label, DM *adaptedDM)
1642: {
1643: PetscBool success;
1647: DMForestTemplate(dm,PetscObjectComm((PetscObject)dm),adaptedDM);
1648: DMForestSetAdaptivityLabel(*adaptedDM,label);
1649: DMSetUp(*adaptedDM);
1650: DMForestGetAdaptivitySuccess(*adaptedDM,&success);
1651: if (!success) {
1652: DMDestroy(adaptedDM);
1653: *adaptedDM = NULL;
1654: }
1655: return(0);
1656: }
1658: static PetscErrorCode DMInitialize_Forest(DM dm)
1659: {
1663: PetscMemzero(dm->ops,sizeof(*(dm->ops)));
1665: dm->ops->clone = DMClone_Forest;
1666: dm->ops->setfromoptions = DMSetFromOptions_Forest;
1667: dm->ops->destroy = DMDestroy_Forest;
1668: dm->ops->createsubdm = DMCreateSubDM_Forest;
1669: dm->ops->refine = DMRefine_Forest;
1670: dm->ops->coarsen = DMCoarsen_Forest;
1671: dm->ops->adaptlabel = DMAdaptLabel_Forest;
1672: return(0);
1673: }
1675: /*MC
1677: DMFOREST = "forest" - A DM object that encapsulates a hierarchically refined mesh. Forests usually have a base DM
1678: (see DMForestGetBaseDM()), from which it is refined. The refinement and partitioning of forests is considered
1679: immutable after DMSetUp() is called. To adapt a mesh, one should call DMForestTemplate() to create a new mesh that
1680: will default to being identical to it, specify how that mesh should differ, and then calling DMSetUp() on the new
1681: mesh.
1683: To specify that a mesh should be refined or coarsened from the previous mesh, a label should be defined on the
1684: previous mesh whose values indicate which cells should be refined (DM_ADAPT_REFINE) or coarsened (DM_ADAPT_COARSEN)
1685: and how (subtypes are free to allow additional values for things like anisotropic refinement). The label should be
1686: given to the *new* mesh with DMForestSetAdaptivityLabel().
1688: Level: advanced
1690: .seealso: DMType, DMCreate(), DMSetType(), DMForestGetBaseDM(), DMForestSetBaseDM(), DMForestTemplate(), DMForestSetAdaptivityLabel()
1691: M*/
1693: PETSC_EXTERN PetscErrorCode DMCreate_Forest(DM dm)
1694: {
1695: DM_Forest *forest;
1700: PetscNewLog(dm,&forest);
1701: dm->dim = 0;
1702: dm->data = forest;
1703: forest->refct = 1;
1704: forest->data = NULL;
1705: forest->topology = NULL;
1706: forest->adapt = NULL;
1707: forest->base = NULL;
1708: forest->adaptPurpose = DM_ADAPT_DETERMINE;
1709: forest->adjDim = PETSC_DEFAULT;
1710: forest->overlap = PETSC_DEFAULT;
1711: forest->minRefinement = PETSC_DEFAULT;
1712: forest->maxRefinement = PETSC_DEFAULT;
1713: forest->initRefinement = PETSC_DEFAULT;
1714: forest->cStart = PETSC_DETERMINE;
1715: forest->cEnd = PETSC_DETERMINE;
1716: forest->cellSF = NULL;
1717: forest->adaptLabel = NULL;
1718: forest->gradeFactor = 2;
1719: forest->cellWeights = NULL;
1720: forest->cellWeightsCopyMode = PETSC_USE_POINTER;
1721: forest->weightsFactor = 1.;
1722: forest->weightCapacity = 1.;
1723: DMForestSetAdaptivityStrategy(dm,DMFORESTADAPTALL);
1724: DMInitialize_Forest(dm);
1725: return(0);
1726: }