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;
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: DMGetMatType(dm,&mtype);
178: DMSetMatType(*tdm,mtype);
179: return(0);
180: }
182: static PetscErrorCode DMInitialize_Forest(DM dm);
184: PETSC_EXTERN PetscErrorCode DMClone_Forest(DM dm, DM *newdm)
185: {
186: DM_Forest *forest = (DM_Forest*) dm->data;
187: const char *type;
191: forest->refct++;
192: (*newdm)->data = forest;
193: PetscObjectGetType((PetscObject) dm, &type);
194: PetscObjectChangeTypeName((PetscObject) *newdm, type);
195: DMInitialize_Forest(*newdm);
196: return(0);
197: }
199: static PetscErrorCode DMDestroy_Forest(DM dm)
200: {
201: DM_Forest *forest = (DM_Forest*) dm->data;
205: if (--forest->refct > 0) return(0);
206: if (forest->destroy) {(*forest->destroy)(dm);}
207: PetscSFDestroy(&forest->cellSF);
208: PetscSFDestroy(&forest->preCoarseToFine);
209: PetscSFDestroy(&forest->coarseToPreFine);
210: DMLabelDestroy(&forest->adaptLabel);
211: PetscFree(forest->adaptStrategy);
212: DMDestroy(&forest->base);
213: DMDestroy(&forest->adapt);
214: PetscFree(forest->topology);
215: PetscFree(forest);
216: return(0);
217: }
219: /*@C
220: DMForestSetTopology - Set the topology of a DMForest during the pre-setup phase. The topology is a string (e.g.
221: "cube", "shell") and can be interpreted by subtypes of DMFOREST) to construct the base DM of a forest during
222: DMSetUp().
224: Logically collective on dm
226: Input parameters:
227: + dm - the forest
228: - topology - the topology of the forest
230: Level: intermediate
232: .seealso: DMForestGetTopology(), DMForestSetBaseDM()
233: @*/
234: PetscErrorCode DMForestSetTopology(DM dm, DMForestTopology topology)
235: {
236: DM_Forest *forest = (DM_Forest*) dm->data;
241: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the topology after setup");
242: PetscFree(forest->topology);
243: PetscStrallocpy((const char*)topology,(char**) &forest->topology);
244: return(0);
245: }
247: /*@C
248: DMForestGetTopology - Get a string describing the topology of a DMForest.
250: Not collective
252: Input parameter:
253: . dm - the forest
255: Output parameter:
256: . topology - the topology of the forest (e.g., 'cube', 'shell')
258: Level: intermediate
260: .seealso: DMForestSetTopology()
261: @*/
262: PetscErrorCode DMForestGetTopology(DM dm, DMForestTopology *topology)
263: {
264: DM_Forest *forest = (DM_Forest*) dm->data;
269: *topology = forest->topology;
270: return(0);
271: }
273: /*@
274: DMForestSetBaseDM - During the pre-setup phase, set the DM that defines the base mesh of a DMForest forest. The
275: forest will be hierarchically refined from the base, and all refinements/coarsenings of the forest will share its
276: base. In general, two forest must share a base to be comparable, to do things like construct interpolators.
278: Logically collective on dm
280: Input Parameters:
281: + dm - the forest
282: - base - the base DM of the forest
284: Notes:
285: Currently the base DM must be a DMPLEX
287: Level: intermediate
289: .seealso: DMForestGetBaseDM()
290: @*/
291: PetscErrorCode DMForestSetBaseDM(DM dm, DM base)
292: {
293: DM_Forest *forest = (DM_Forest*) dm->data;
294: PetscInt dim, dimEmbed;
299: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the base after setup");
300: PetscObjectReference((PetscObject)base);
301: DMDestroy(&forest->base);
302: forest->base = base;
303: if (base) {
304: PetscBool isper;
305: const PetscReal *maxCell, *L;
306: const DMBoundaryType *bd;
309: DMGetDimension(base,&dim);
310: DMSetDimension(dm,dim);
311: DMGetCoordinateDim(base,&dimEmbed);
312: DMSetCoordinateDim(dm,dimEmbed);
313: DMGetPeriodicity(base,&isper,&maxCell,&L,&bd);
314: DMSetPeriodicity(dm,isper,maxCell,L,bd);
315: } else {
316: DMSetPeriodicity(dm,PETSC_FALSE,NULL,NULL,NULL);
317: }
318: return(0);
319: }
321: /*@
322: DMForestGetBaseDM - Get the base DM of a DMForest forest. The forest will be hierarchically refined from the base,
323: and all refinements/coarsenings of the forest will share its base. In general, two forest must share a base to be
324: comparable, to do things like construct interpolators.
326: Not collective
328: Input Parameter:
329: . dm - the forest
331: Output Parameter:
332: . base - the base DM of the forest
334: Notes:
335: After DMSetUp(), the base DM will be redundantly distributed across MPI processes
337: Level: intermediate
339: .seealso: DMForestSetBaseDM()
340: @*/
341: PetscErrorCode DMForestGetBaseDM(DM dm, DM *base)
342: {
343: DM_Forest *forest = (DM_Forest*) dm->data;
348: *base = forest->base;
349: return(0);
350: }
352: PetscErrorCode DMForestSetBaseCoordinateMapping(DM dm, PetscErrorCode (*func)(DM,PetscInt,PetscInt,const PetscReal [],PetscReal [],void*),void *ctx)
353: {
354: DM_Forest *forest = (DM_Forest*) dm->data;
358: forest->mapcoordinates = func;
359: forest->mapcoordinatesctx = ctx;
360: return(0);
361: }
363: PetscErrorCode DMForestGetBaseCoordinateMapping(DM dm, PetscErrorCode (**func) (DM,PetscInt,PetscInt,const PetscReal [],PetscReal [],void*),void *ctx)
364: {
365: DM_Forest *forest = (DM_Forest*) dm->data;
369: if (func) *func = forest->mapcoordinates;
370: if (ctx) *((void**) ctx) = forest->mapcoordinatesctx;
371: return(0);
372: }
374: /*@
375: DMForestSetAdaptivityForest - During the pre-setup phase, set the forest from which the current forest will be
376: adapted (e.g., the current forest will be refined/coarsened/repartitioned from it) im DMSetUp(). Usually not needed
377: by users directly: DMForestTemplate() constructs a new forest to be adapted from an old forest and calls this
378: routine.
380: Note that this can be called after setup with adapt = NULL, which will clear all internal data related to the
381: adaptivity forest from dm. This way, repeatedly adapting does not leave stale DM objects in memory.
383: Logically collective on dm
385: Input Parameters:
386: + dm - the new forest, which will be constructed from adapt
387: - adapt - the old forest
389: Level: intermediate
391: .seealso: DMForestGetAdaptivityForest(), DMForestSetAdaptivityPurpose()
392: @*/
393: PetscErrorCode DMForestSetAdaptivityForest(DM dm,DM adapt)
394: {
395: DM_Forest *forest, *adaptForest, *oldAdaptForest;
396: DM oldAdapt;
397: PetscBool isForest;
403: DMIsForest(dm, &isForest);
404: if (!isForest) return(0);
405: if (adapt != NULL && dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the adaptation forest after setup");
406: forest = (DM_Forest*) dm->data;
407: DMForestGetAdaptivityForest(dm,&oldAdapt);
408: adaptForest = (DM_Forest*) (adapt ? adapt->data : NULL);
409: oldAdaptForest = (DM_Forest*) (oldAdapt ? oldAdapt->data : NULL);
410: if (adaptForest != oldAdaptForest) {
411: PetscSFDestroy(&forest->preCoarseToFine);
412: PetscSFDestroy(&forest->coarseToPreFine);
413: if (forest->clearadaptivityforest) {(*forest->clearadaptivityforest)(dm);}
414: }
415: switch (forest->adaptPurpose) {
416: case DM_ADAPT_DETERMINE:
417: PetscObjectReference((PetscObject)adapt);
418: DMDestroy(&(forest->adapt));
419: forest->adapt = adapt;
420: break;
421: case DM_ADAPT_REFINE:
422: DMSetCoarseDM(dm,adapt);
423: break;
424: case DM_ADAPT_COARSEN:
425: case DM_ADAPT_COARSEN_LAST:
426: DMSetFineDM(dm,adapt);
427: break;
428: default:
429: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"invalid adaptivity purpose");
430: }
431: return(0);
432: }
434: /*@
435: DMForestGetAdaptivityForest - Get the forest from which the current forest is adapted.
437: Not collective
439: Input Parameter:
440: . dm - the forest
442: Output Parameter:
443: . adapt - the forest from which dm is/was adapted
445: Level: intermediate
447: .seealso: DMForestSetAdaptivityForest(), DMForestSetAdaptivityPurpose()
448: @*/
449: PetscErrorCode DMForestGetAdaptivityForest(DM dm, DM *adapt)
450: {
451: DM_Forest *forest;
456: forest = (DM_Forest*) dm->data;
457: switch (forest->adaptPurpose) {
458: case DM_ADAPT_DETERMINE:
459: *adapt = forest->adapt;
460: break;
461: case DM_ADAPT_REFINE:
462: DMGetCoarseDM(dm,adapt);
463: break;
464: case DM_ADAPT_COARSEN:
465: case DM_ADAPT_COARSEN_LAST:
466: DMGetFineDM(dm,adapt);
467: break;
468: default:
469: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"invalid adaptivity purpose");
470: }
471: return(0);
472: }
474: /*@
475: DMForestSetAdaptivityPurpose - During the pre-setup phase, set whether the current DM is being adapted from its
476: source (set with DMForestSetAdaptivityForest()) for the purpose of refinement (DM_ADAPT_REFINE), coarsening
477: (DM_ADAPT_COARSEN), or undefined (DM_ADAPT_DETERMINE). This only matters for the purposes of reference counting:
478: during DMDestroy(), cyclic references can be found between DMs only if the cyclic reference is due to a fine/coarse
479: relationship (see DMSetFineDM()/DMSetCoarseDM()). If the purpose is not refinement or coarsening, and the user does
480: not maintain a reference to the post-adaptation forest (i.e., the one created by DMForestTemplate()), then this can
481: cause a memory leak. This method is used by subtypes of DMForest when automatically constructing mesh hierarchies.
483: Logically collective on dm
485: Input Parameters:
486: + dm - the forest
487: - purpose - the adaptivity purpose
489: Level: advanced
491: .seealso: DMForestTemplate(), DMForestSetAdaptivityForest(), DMForestGetAdaptivityForest(), DMAdaptFlag
492: @*/
493: PetscErrorCode DMForestSetAdaptivityPurpose(DM dm, DMAdaptFlag purpose)
494: {
495: DM_Forest *forest;
499: forest = (DM_Forest*) dm->data;
500: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the adaptation forest after setup");
501: if (purpose != forest->adaptPurpose) {
502: DM adapt;
504: DMForestGetAdaptivityForest(dm,&adapt);
505: PetscObjectReference((PetscObject)adapt);
506: DMForestSetAdaptivityForest(dm,NULL);
508: forest->adaptPurpose = purpose;
510: DMForestSetAdaptivityForest(dm,adapt);
511: DMDestroy(&adapt);
512: }
513: return(0);
514: }
516: /*@
517: DMForestGetAdaptivityPurpose - Get whether the current DM is being adapted from its source (set with
518: DMForestSetAdaptivityForest()) for the purpose of refinement (DM_ADAPT_REFINE), coarsening (DM_ADAPT_COARSEN),
519: coarsening only the last level (DM_ADAPT_COARSEN_LAST) or undefined (DM_ADAPT_DETERMINE).
520: This only matters for the purposes of reference counting: during DMDestroy(), cyclic
521: references can be found between DMs only if the cyclic reference is due to a fine/coarse relationship (see
522: DMSetFineDM()/DMSetCoarseDM()). If the purpose is not refinement or coarsening, and the user does not maintain a
523: reference to the post-adaptation forest (i.e., the one created by DMForestTemplate()), then this can cause a memory
524: leak. This method is used by subtypes of DMForest when automatically constructing mesh hierarchies.
526: Not collective
528: Input Parameter:
529: . dm - the forest
531: Output Parameter:
532: . purpose - the adaptivity purpose
534: Level: advanced
536: .seealso: DMForestTemplate(), DMForestSetAdaptivityForest(), DMForestGetAdaptivityForest(), DMAdaptFlag
537: @*/
538: PetscErrorCode DMForestGetAdaptivityPurpose(DM dm, DMAdaptFlag *purpose)
539: {
540: DM_Forest *forest;
543: forest = (DM_Forest*) dm->data;
544: *purpose = forest->adaptPurpose;
545: return(0);
546: }
548: /*@
549: DMForestSetAdjacencyDimension - During the pre-setup phase, set the dimension of interface points that determine
550: cell adjacency (for the purposes of partitioning and overlap).
552: Logically collective on dm
554: Input Parameters:
555: + dm - the forest
556: - adjDim - default 0 (i.e., vertices determine adjacency)
558: Level: intermediate
560: .seealso: DMForestGetAdjacencyDimension(), DMForestSetAdjacencyCodimension(), DMForestSetPartitionOverlap()
561: @*/
562: PetscErrorCode DMForestSetAdjacencyDimension(DM dm, PetscInt adjDim)
563: {
564: PetscInt dim;
565: DM_Forest *forest = (DM_Forest*) dm->data;
570: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the adjacency dimension after setup");
571: if (adjDim < 0) SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_OUTOFRANGE,"adjacency dim cannot be < 0: %d", adjDim);
572: DMGetDimension(dm,&dim);
573: if (adjDim > dim) SETERRQ2(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_OUTOFRANGE,"adjacency dim cannot be > %d: %d", dim, adjDim);
574: forest->adjDim = adjDim;
575: return(0);
576: }
578: /*@
579: DMForestSetAdjacencyCodimension - Like DMForestSetAdjacencyDimension(), but specified as a co-dimension (so that,
580: e.g., adjacency based on facets can be specified by codimension 1 in all cases)
582: Logically collective on dm
584: Input Parameters:
585: + dm - the forest
586: - adjCodim - default isthe dimension of the forest (see DMGetDimension()), since this is the codimension of vertices
588: Level: intermediate
590: .seealso: DMForestGetAdjacencyCodimension(), DMForestSetAdjacencyDimension()
591: @*/
592: PetscErrorCode DMForestSetAdjacencyCodimension(DM dm, PetscInt adjCodim)
593: {
594: PetscInt dim;
599: DMGetDimension(dm,&dim);
600: DMForestSetAdjacencyDimension(dm,dim-adjCodim);
601: return(0);
602: }
604: /*@
605: DMForestGetAdjacencyDimension - Get the dimension of interface points that determine cell adjacency (for the
606: purposes of partitioning and overlap).
608: Not collective
610: Input Parameter:
611: . dm - the forest
613: Output Parameter:
614: . adjDim - default 0 (i.e., vertices determine adjacency)
616: Level: intermediate
618: .seealso: DMForestSetAdjacencyDimension(), DMForestGetAdjacencyCodimension(), DMForestSetPartitionOverlap()
619: @*/
620: PetscErrorCode DMForestGetAdjacencyDimension(DM dm, PetscInt *adjDim)
621: {
622: DM_Forest *forest = (DM_Forest*) dm->data;
627: *adjDim = forest->adjDim;
628: return(0);
629: }
631: /*@
632: DMForestGetAdjacencyCodimension - Like DMForestGetAdjacencyDimension(), but specified as a co-dimension (so that,
633: e.g., adjacency based on facets can be specified by codimension 1 in all cases)
635: Not collective
637: Input Parameter:
638: . dm - the forest
640: Output Parameter:
641: . adjCodim - default isthe dimension of the forest (see DMGetDimension()), since this is the codimension of vertices
643: Level: intermediate
645: .seealso: DMForestSetAdjacencyCodimension(), DMForestGetAdjacencyDimension()
646: @*/
647: PetscErrorCode DMForestGetAdjacencyCodimension(DM dm, PetscInt *adjCodim)
648: {
649: DM_Forest *forest = (DM_Forest*) dm->data;
650: PetscInt dim;
656: DMGetDimension(dm,&dim);
657: *adjCodim = dim - forest->adjDim;
658: return(0);
659: }
661: /*@
662: DMForestSetPartitionOverlap - During the pre-setup phase, set the amount of cell-overlap present in parallel
663: partitions of a forest, with values > 0 indicating subdomains that are expanded by that many iterations of adding
664: adjacent cells
666: Logically collective on dm
668: Input Parameters:
669: + dm - the forest
670: - overlap - default 0
672: Level: intermediate
674: .seealso: DMForestGetPartitionOverlap(), DMForestSetAdjacencyDimension(), DMForestSetAdjacencyCodimension()
675: @*/
676: PetscErrorCode DMForestSetPartitionOverlap(DM dm, PetscInt overlap)
677: {
678: DM_Forest *forest = (DM_Forest*) dm->data;
682: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the overlap after setup");
683: if (overlap < 0) SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_OUTOFRANGE,"overlap cannot be < 0: %d", overlap);
684: forest->overlap = overlap;
685: return(0);
686: }
688: /*@
689: DMForestGetPartitionOverlap - Get the amount of cell-overlap present in parallel partitions of a forest, with values
690: > 0 indicating subdomains that are expanded by that many iterations of adding adjacent cells
692: Not collective
694: Input Parameter:
695: . dm - the forest
697: Output Parameter:
698: . overlap - default 0
700: Level: intermediate
702: .seealso: DMForestGetPartitionOverlap(), DMForestSetAdjacencyDimension(), DMForestSetAdjacencyCodimension()
703: @*/
704: PetscErrorCode DMForestGetPartitionOverlap(DM dm, PetscInt *overlap)
705: {
706: DM_Forest *forest = (DM_Forest*) dm->data;
711: *overlap = forest->overlap;
712: return(0);
713: }
715: /*@
716: DMForestSetMinimumRefinement - During the pre-setup phase, set the minimum level of refinement (relative to the base
717: DM, see DMForestGetBaseDM()) allowed in the forest. If the forest is being created by coarsening a previous forest
718: (see DMForestGetAdaptivityForest()) this limits the amount of coarsening.
720: Logically collective on dm
722: Input Parameters:
723: + dm - the forest
724: - minRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
726: Level: intermediate
728: .seealso: DMForestGetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestSetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityForest()
729: @*/
730: PetscErrorCode DMForestSetMinimumRefinement(DM dm, PetscInt minRefinement)
731: {
732: DM_Forest *forest = (DM_Forest*) dm->data;
736: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the minimum refinement after setup");
737: forest->minRefinement = minRefinement;
738: return(0);
739: }
741: /*@
742: DMForestGetMinimumRefinement - Get the minimum level of refinement (relative to the base DM, see
743: DMForestGetBaseDM()) allowed in the forest. If the forest is being created by coarsening a previous forest (see
744: DMForestGetAdaptivityForest()), this limits the amount of coarsening.
746: Not collective
748: Input Parameter:
749: . dm - the forest
751: Output Parameter:
752: . minRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
754: Level: intermediate
756: .seealso: DMForestSetMinimumRefinement(), DMForestGetMaximumRefinement(), DMForestGetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityForest()
757: @*/
758: PetscErrorCode DMForestGetMinimumRefinement(DM dm, PetscInt *minRefinement)
759: {
760: DM_Forest *forest = (DM_Forest*) dm->data;
765: *minRefinement = forest->minRefinement;
766: return(0);
767: }
769: /*@
770: DMForestSetInitialRefinement - During the pre-setup phase, set the initial level of refinement (relative to the base
771: DM, see DMForestGetBaseDM()) allowed in the forest.
773: Logically collective on dm
775: Input Parameters:
776: + dm - the forest
777: - initefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
779: Level: intermediate
781: .seealso: DMForestSetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestGetBaseDM()
782: @*/
783: PetscErrorCode DMForestSetInitialRefinement(DM dm, PetscInt initRefinement)
784: {
785: DM_Forest *forest = (DM_Forest*) dm->data;
789: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the initial refinement after setup");
790: forest->initRefinement = initRefinement;
791: return(0);
792: }
794: /*@
795: DMForestGetInitialRefinement - Get the initial level of refinement (relative to the base DM, see
796: DMForestGetBaseDM()) allowed in the forest.
798: Not collective
800: Input Parameter:
801: . dm - the forest
803: Output Parameter:
804: . initRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
806: Level: intermediate
808: .seealso: DMForestSetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestGetBaseDM()
809: @*/
810: PetscErrorCode DMForestGetInitialRefinement(DM dm, PetscInt *initRefinement)
811: {
812: DM_Forest *forest = (DM_Forest*) dm->data;
817: *initRefinement = forest->initRefinement;
818: return(0);
819: }
821: /*@
822: DMForestSetMaximumRefinement - During the pre-setup phase, set the maximum level of refinement (relative to the base
823: DM, see DMForestGetBaseDM()) allowed in the forest. If the forest is being created by refining a previous forest
824: (see DMForestGetAdaptivityForest()), this limits the amount of refinement.
826: Logically collective on dm
828: Input Parameters:
829: + dm - the forest
830: - maxRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
832: Level: intermediate
834: .seealso: DMForestGetMinimumRefinement(), DMForestSetMaximumRefinement(), DMForestSetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityDM()
835: @*/
836: PetscErrorCode DMForestSetMaximumRefinement(DM dm, PetscInt maxRefinement)
837: {
838: DM_Forest *forest = (DM_Forest*) dm->data;
842: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the maximum refinement after setup");
843: forest->maxRefinement = maxRefinement;
844: return(0);
845: }
847: /*@
848: DMForestGetMaximumRefinement - Get the maximum level of refinement (relative to the base DM, see
849: DMForestGetBaseDM()) allowed in the forest. If the forest is being created by refining a previous forest (see
850: DMForestGetAdaptivityForest()), this limits the amount of refinement.
852: Not collective
854: Input Parameter:
855: . dm - the forest
857: Output Parameter:
858: . maxRefinement - default PETSC_DEFAULT (interpreted by the subtype of DMForest)
860: Level: intermediate
862: .seealso: DMForestSetMaximumRefinement(), DMForestGetMinimumRefinement(), DMForestGetInitialRefinement(), DMForestGetBaseDM(), DMForestGetAdaptivityForest()
863: @*/
864: PetscErrorCode DMForestGetMaximumRefinement(DM dm, PetscInt *maxRefinement)
865: {
866: DM_Forest *forest = (DM_Forest*) dm->data;
871: *maxRefinement = forest->maxRefinement;
872: return(0);
873: }
875: /*@C
876: DMForestSetAdaptivityStrategy - During the pre-setup phase, set the strategy for combining adaptivity labels from multiple processes.
877: Subtypes of DMForest may define their own strategies. Two default strategies are DMFORESTADAPTALL, which indicates that all processes must agree
878: for a refinement/coarsening flag to be valid, and DMFORESTADAPTANY, which indicates that only one process needs to
879: specify refinement/coarsening.
881: Logically collective on dm
883: Input Parameters:
884: + dm - the forest
885: - adaptStrategy - default DMFORESTADAPTALL
887: Level: advanced
889: .seealso: DMForestGetAdaptivityStrategy()
890: @*/
891: PetscErrorCode DMForestSetAdaptivityStrategy(DM dm, DMForestAdaptivityStrategy adaptStrategy)
892: {
893: DM_Forest *forest = (DM_Forest*) dm->data;
898: PetscFree(forest->adaptStrategy);
899: PetscStrallocpy((const char*) adaptStrategy,(char**)&forest->adaptStrategy);
900: return(0);
901: }
903: /*@C
904: DMForestSetAdaptivityStrategy - Get the strategy for combining adaptivity labels from multiple processes. Subtypes
905: of DMForest may define their own strategies. Two default strategies are DMFORESTADAPTALL, which indicates that all
906: processes must agree for a refinement/coarsening flag to be valid, and DMFORESTADAPTANY, which indicates that only
907: one process needs to specify refinement/coarsening.
909: Not collective
911: Input Parameter:
912: . dm - the forest
914: Output Parameter:
915: . adaptStrategy - the adaptivity strategy (default DMFORESTADAPTALL)
917: Level: advanced
919: .seealso: DMForestSetAdaptivityStrategy()
920: @*/
921: PetscErrorCode DMForestGetAdaptivityStrategy(DM dm, DMForestAdaptivityStrategy *adaptStrategy)
922: {
923: DM_Forest *forest = (DM_Forest*) dm->data;
928: *adaptStrategy = forest->adaptStrategy;
929: return(0);
930: }
932: /*@
933: DMForestGetAdaptivitySuccess - Return whether the requested adaptation (refinement, coarsening, repartitioning,
934: etc.) was successful. PETSC_FALSE indicates that the post-adaptation forest is the same as the pre-adpatation
935: forest. A requested adaptation may have been unsuccessful if, for example, the requested refinement would have
936: exceeded the maximum refinement level.
938: Collective on dm
940: Input Parameter:
942: . dm - the post-adaptation forest
944: Output Parameter:
946: . success - PETSC_TRUE if the post-adaptation forest is different from the pre-adaptation forest.
948: Level: intermediate
950: .see
951: @*/
952: PetscErrorCode DMForestGetAdaptivitySuccess(DM dm, PetscBool *success)
953: {
954: DM_Forest *forest;
959: if (!dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"DMSetUp() has not been called yet.");
960: forest = (DM_Forest *) dm->data;
961: (forest->getadaptivitysuccess)(dm,success);
962: return(0);
963: }
965: /*@
966: DMForestSetComputeAdaptivitySF - During the pre-setup phase, set whether transfer PetscSFs should be computed
967: 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().
969: Logically collective on dm
971: Input Parameters:
972: + dm - the post-adaptation forest
973: - computeSF - default PETSC_TRUE
975: Level: advanced
977: .seealso: DMForestGetComputeAdaptivitySF(), DMForestGetAdaptivitySF()
978: @*/
979: PetscErrorCode DMForestSetComputeAdaptivitySF(DM dm, PetscBool computeSF)
980: {
981: DM_Forest *forest;
985: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot compute adaptivity PetscSFs after setup is called");
986: forest = (DM_Forest*) dm->data;
987: forest->computeAdaptSF = computeSF;
988: return(0);
989: }
991: PetscErrorCode DMForestTransferVec(DM dmIn, Vec vecIn, DM dmOut, Vec vecOut, PetscBool useBCs, PetscReal time)
992: {
993: DM_Forest *forest;
1001: forest = (DM_Forest *) dmIn->data;
1002: if (!forest->transfervec) SETERRQ(PetscObjectComm((PetscObject)dmIn),PETSC_ERR_SUP,"DMForestTransferVec() not implemented");
1003: (forest->transfervec)(dmIn,vecIn,dmOut,vecOut,useBCs,time);
1004: return(0);
1005: }
1007: PetscErrorCode DMForestTransferVecFromBase(DM dm, Vec vecIn, Vec vecOut)
1008: {
1009: DM_Forest *forest;
1016: forest = (DM_Forest *) dm->data;
1017: if (!forest->transfervecfrombase) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"DMForestTransferVecFromBase() not implemented");
1018: (forest->transfervecfrombase)(dm,vecIn,vecOut);
1019: return(0);
1020: }
1022: /*@
1023: DMForestGetComputeAdaptivitySF - Get whether transfer PetscSFs should be computed relating the cells of the
1024: pre-adaptation forest to the post-adaptiation forest. After DMSetUp() is called, these transfer PetscSFs can be
1025: accessed with DMForestGetAdaptivitySF().
1027: Not collective
1029: Input Parameter:
1030: . dm - the post-adaptation forest
1032: Output Parameter:
1033: . computeSF - default PETSC_TRUE
1035: Level: advanced
1037: .seealso: DMForestSetComputeAdaptivitySF(), DMForestGetAdaptivitySF()
1038: @*/
1039: PetscErrorCode DMForestGetComputeAdaptivitySF(DM dm, PetscBool *computeSF)
1040: {
1041: DM_Forest *forest;
1045: forest = (DM_Forest*) dm->data;
1046: *computeSF = forest->computeAdaptSF;
1047: return(0);
1048: }
1050: /*@
1051: DMForestGetAdaptivitySF - Get PetscSFs that relate the pre-adaptation forest to the post-adaptation forest.
1052: Adaptation can be any combination of refinement, coarsening, repartition, and change of overlap, so there may be
1053: some cells of the pre-adaptation that are parents of post-adaptation cells, and vice versa. Therefore there are two
1054: PetscSFs: one that relates pre-adaptation coarse cells to post-adaptation fine cells, and one that relates
1055: pre-adaptation fine cells to post-adaptation coarse cells.
1057: Not collective
1059: Input Parameter:
1060: dm - the post-adaptation forest
1062: Output Parameter:
1063: preCoarseToFine - pre-adaptation coarse cells to post-adaptation fine cells: BCast goes from pre- to post-
1064: coarseToPreFine - post-adaptation coarse cells to pre-adaptation fine cells: BCast goes from post- to pre-
1066: Level: advanced
1068: .seealso: DMForestGetComputeAdaptivitySF(), DMForestSetComputeAdaptivitySF()
1069: @*/
1070: PetscErrorCode DMForestGetAdaptivitySF(DM dm, PetscSF *preCoarseToFine, PetscSF *coarseToPreFine)
1071: {
1072: DM_Forest *forest;
1077: DMSetUp(dm);
1078: forest = (DM_Forest*) dm->data;
1079: if (preCoarseToFine) *preCoarseToFine = forest->preCoarseToFine;
1080: if (coarseToPreFine) *coarseToPreFine = forest->coarseToPreFine;
1081: return(0);
1082: }
1084: /*@
1085: DMForestSetGradeFactor - During the pre-setup phase, set the desired amount of grading in the mesh, e.g. give 2 to
1086: indicate that the diameter of neighboring cells should differ by at most a factor of 2. Subtypes of DMForest may
1087: only support one particular choice of grading factor.
1089: Logically collective on dm
1091: Input Parameters:
1092: + dm - the forest
1093: - grade - the grading factor
1095: Level: advanced
1097: .seealso: DMForestGetGradeFactor()
1098: @*/
1099: PetscErrorCode DMForestSetGradeFactor(DM dm, PetscInt grade)
1100: {
1101: DM_Forest *forest = (DM_Forest*) dm->data;
1105: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the grade factor after setup");
1106: forest->gradeFactor = grade;
1107: return(0);
1108: }
1110: /*@
1111: DMForestGetGradeFactor - Get the desired amount of grading in the mesh, e.g. give 2 to indicate that the diameter of
1112: neighboring cells should differ by at most a factor of 2. Subtypes of DMForest may only support one particular
1113: choice of grading factor.
1115: Not collective
1117: Input Parameter:
1118: . dm - the forest
1120: Output Parameter:
1121: . grade - the grading factor
1123: Level: advanced
1125: .seealso: DMForestSetGradeFactor()
1126: @*/
1127: PetscErrorCode DMForestGetGradeFactor(DM dm, PetscInt *grade)
1128: {
1129: DM_Forest *forest = (DM_Forest*) dm->data;
1134: *grade = forest->gradeFactor;
1135: return(0);
1136: }
1138: /*@
1139: DMForestSetCellWeightFactor - During the pre-setup phase, set the factor by which the level of refinement changes
1140: the cell weight (see DMForestSetCellWeights()) when calculating partitions. The final weight of a cell will be
1141: (cellWeight) * (weightFactor^refinementLevel). A factor of 1 indicates that the weight of a cell does not depend on
1142: its level; a factor of 2, for example, might be appropriate for sub-cycling time-stepping methods, when the
1143: computation associated with a cell is multiplied by a factor of 2 for each additional level of refinement.
1145: Logically collective on dm
1147: Input Parameters:
1148: + dm - the forest
1149: - weightsFactors - default 1.
1151: Level: advanced
1153: .seealso: DMForestGetCellWeightFactor(), DMForestSetCellWeights()
1154: @*/
1155: PetscErrorCode DMForestSetCellWeightFactor(DM dm, PetscReal weightsFactor)
1156: {
1157: DM_Forest *forest = (DM_Forest*) dm->data;
1161: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the weights factor after setup");
1162: forest->weightsFactor = weightsFactor;
1163: return(0);
1164: }
1166: /*@
1167: DMForestGetCellWeightFactor - Get the factor by which the level of refinement changes the cell weight (see
1168: DMForestSetCellWeights()) when calculating partitions. The final weight of a cell will be (cellWeight) *
1169: (weightFactor^refinementLevel). A factor of 1 indicates that the weight of a cell does not depend on its level; a
1170: factor of 2, for example, might be appropriate for sub-cycling time-stepping methods, when the computation
1171: associated with a cell is multiplied by a factor of 2 for each additional level of refinement.
1173: Not collective
1175: Input Parameter:
1176: . dm - the forest
1178: Output Parameter:
1179: . weightsFactors - default 1.
1181: Level: advanced
1183: .seealso: DMForestSetCellWeightFactor(), DMForestSetCellWeights()
1184: @*/
1185: PetscErrorCode DMForestGetCellWeightFactor(DM dm, PetscReal *weightsFactor)
1186: {
1187: DM_Forest *forest = (DM_Forest*) dm->data;
1192: *weightsFactor = forest->weightsFactor;
1193: return(0);
1194: }
1196: /*@
1197: DMForestGetCellChart - After the setup phase, get the local half-open interval of the chart of cells on this process
1199: Not collective
1201: Input Parameter:
1202: . dm - the forest
1204: Output Parameters:
1205: + cStart - the first cell on this process
1206: - cEnd - one after the final cell on this process
1208: Level: intermediate
1210: .seealso: DMForestGetCellSF()
1211: @*/
1212: PetscErrorCode DMForestGetCellChart(DM dm, PetscInt *cStart, PetscInt *cEnd)
1213: {
1214: DM_Forest *forest = (DM_Forest*) dm->data;
1221: if (((forest->cStart == PETSC_DETERMINE) || (forest->cEnd == PETSC_DETERMINE)) && forest->createcellchart) {
1222: forest->createcellchart(dm,&forest->cStart,&forest->cEnd);
1223: }
1224: *cStart = forest->cStart;
1225: *cEnd = forest->cEnd;
1226: return(0);
1227: }
1229: /*@
1230: DMForestGetCellSF - After the setup phase, get the PetscSF for overlapping cells between processes
1232: Not collective
1234: Input Parameter:
1235: . dm - the forest
1237: Output Parameter:
1238: . cellSF - the PetscSF
1240: Level: intermediate
1242: .seealso: DMForestGetCellChart()
1243: @*/
1244: PetscErrorCode DMForestGetCellSF(DM dm, PetscSF *cellSF)
1245: {
1246: DM_Forest *forest = (DM_Forest*) dm->data;
1252: if ((!forest->cellSF) && forest->createcellsf) {
1253: forest->createcellsf(dm,&forest->cellSF);
1254: }
1255: *cellSF = forest->cellSF;
1256: return(0);
1257: }
1259: /*@C
1260: DMForestSetAdaptivityLabel - During the pre-setup phase, set the label of the pre-adaptation forest (see
1261: DMForestGetAdaptivityForest()) that holds the adaptation flags (refinement, coarsening, or some combination). The
1262: interpretation of the label values is up to the subtype of DMForest, but DM_ADAPT_DETERMINE, DM_ADAPT_KEEP,
1263: DM_ADAPT_REFINE, and DM_ADAPT_COARSEN have been reserved as choices that should be accepted by all subtypes.
1265: Logically collective on dm
1267: Input Parameters:
1268: - dm - the forest
1269: + adaptLabel - the label in the pre-adaptation forest
1271: Level: intermediate
1273: .seealso DMForestGetAdaptivityLabel()
1274: @*/
1275: PetscErrorCode DMForestSetAdaptivityLabel(DM dm, DMLabel adaptLabel)
1276: {
1277: DM_Forest *forest = (DM_Forest*) dm->data;
1283: PetscObjectReference((PetscObject)adaptLabel);
1284: DMLabelDestroy(&forest->adaptLabel);
1285: forest->adaptLabel = adaptLabel;
1286: return(0);
1287: }
1289: /*@C
1290: DMForestGetAdaptivityLabel - Get the label of the pre-adaptation forest (see DMForestGetAdaptivityForest()) that
1291: holds the adaptation flags (refinement, coarsening, or some combination). The interpretation of the label values is
1292: up to the subtype of DMForest, but DM_ADAPT_DETERMINE, DM_ADAPT_KEEP, DM_ADAPT_REFINE, and DM_ADAPT_COARSEN have
1293: been reserved as choices that should be accepted by all subtypes.
1295: Not collective
1297: Input Parameter:
1298: . dm - the forest
1300: Output Parameter:
1301: . adaptLabel - the name of the label in the pre-adaptation forest
1303: Level: intermediate
1305: .seealso DMForestSetAdaptivityLabel()
1306: @*/
1307: PetscErrorCode DMForestGetAdaptivityLabel(DM dm, DMLabel *adaptLabel)
1308: {
1309: DM_Forest *forest = (DM_Forest*) dm->data;
1313: *adaptLabel = forest->adaptLabel;
1314: return(0);
1315: }
1317: /*@
1318: DMForestSetCellWeights - Set the weights assigned to each of the cells (see DMForestGetCellChart()) of the current
1319: process: weights are used to determine parallel partitioning. Partitions will be created so that each process's
1320: ratio of weight to capacity (see DMForestSetWeightCapacity()) is roughly equal. If NULL, each cell receives a weight
1321: of 1.
1323: Logically collective on dm
1325: Input Parameters:
1326: + dm - the forest
1327: . weights - the array of weights for all cells, or NULL to indicate each cell has weight 1.
1328: - copyMode - how weights should reference weights
1330: Level: advanced
1332: .seealso: DMForestGetCellWeights(), DMForestSetWeightCapacity()
1333: @*/
1334: PetscErrorCode DMForestSetCellWeights(DM dm, PetscReal weights[], PetscCopyMode copyMode)
1335: {
1336: DM_Forest *forest = (DM_Forest*) dm->data;
1337: PetscInt cStart, cEnd;
1342: DMForestGetCellChart(dm,&cStart,&cEnd);
1343: if (cEnd < cStart) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"cell chart [%d,%d) is not valid",cStart,cEnd);
1344: if (copyMode == PETSC_COPY_VALUES) {
1345: if (forest->cellWeightsCopyMode != PETSC_OWN_POINTER || forest->cellWeights == weights) {
1346: PetscMalloc1(cEnd-cStart,&forest->cellWeights);
1347: }
1348: PetscArraycpy(forest->cellWeights,weights,cEnd-cStart);
1349: forest->cellWeightsCopyMode = PETSC_OWN_POINTER;
1350: return(0);
1351: }
1352: if (forest->cellWeightsCopyMode == PETSC_OWN_POINTER) {
1353: PetscFree(forest->cellWeights);
1354: }
1355: forest->cellWeights = weights;
1356: forest->cellWeightsCopyMode = copyMode;
1357: return(0);
1358: }
1360: /*@
1361: DMForestGetCellWeights - Get the weights assigned to each of the cells (see DMForestGetCellChart()) of the current
1362: process: weights are used to determine parallel partitioning. Partitions will be created so that each process's
1363: ratio of weight to capacity (see DMForestSetWeightCapacity()) is roughly equal. If NULL, each cell receives a weight
1364: of 1.
1366: Not collective
1368: Input Parameter:
1369: . dm - the forest
1371: Output Parameter:
1372: . weights - the array of weights for all cells, or NULL to indicate each cell has weight 1.
1374: Level: advanced
1376: .seealso: DMForestSetCellWeights(), DMForestSetWeightCapacity()
1377: @*/
1378: PetscErrorCode DMForestGetCellWeights(DM dm, PetscReal **weights)
1379: {
1380: DM_Forest *forest = (DM_Forest*) dm->data;
1385: *weights = forest->cellWeights;
1386: return(0);
1387: }
1389: /*@
1390: DMForestSetWeightCapacity - During the pre-setup phase, set the capacity of the current process when repartitioning
1391: a pre-adaptation forest (see DMForestGetAdaptivityForest()). After partitioning, the ratio of the weight of each
1392: process's cells to the process's capacity will be roughly equal for all processes. A capacity of 0 indicates that
1393: the current process should not have any cells after repartitioning.
1395: Logically Collective on dm
1397: Input parameters:
1398: + dm - the forest
1399: - capacity - this process's capacity
1401: Level: advanced
1403: .seealso DMForestGetWeightCapacity(), DMForestSetCellWeights(), DMForestSetCellWeightFactor()
1404: @*/
1405: PetscErrorCode DMForestSetWeightCapacity(DM dm, PetscReal capacity)
1406: {
1407: DM_Forest *forest = (DM_Forest*) dm->data;
1411: if (dm->setupcalled) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_WRONGSTATE,"Cannot change the weight capacity after setup");
1412: if (capacity < 0.) SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_OUTOFRANGE,"Cannot have negative weight capacity; %f",capacity);
1413: forest->weightCapacity = capacity;
1414: return(0);
1415: }
1417: /*@
1418: DMForestGetWeightCapacity - Set the capacity of the current process when repartitioning a pre-adaptation forest (see
1419: DMForestGetAdaptivityForest()). After partitioning, the ratio of the weight of each process's cells to the
1420: process's capacity will be roughly equal for all processes. A capacity of 0 indicates that the current process
1421: should not have any cells after repartitioning.
1423: Not collective
1425: Input parameter:
1426: . dm - the forest
1428: Output parameter:
1429: . capacity - this process's capacity
1431: Level: advanced
1433: .seealso DMForestSetWeightCapacity(), DMForestSetCellWeights(), DMForestSetCellWeightFactor()
1434: @*/
1435: PetscErrorCode DMForestGetWeightCapacity(DM dm, PetscReal *capacity)
1436: {
1437: DM_Forest *forest = (DM_Forest*) dm->data;
1442: *capacity = forest->weightCapacity;
1443: return(0);
1444: }
1446: PETSC_EXTERN PetscErrorCode DMSetFromOptions_Forest(PetscOptionItems *PetscOptionsObject,DM dm)
1447: {
1448: PetscBool flg, flg1, flg2, flg3, flg4;
1449: DMForestTopology oldTopo;
1450: char stringBuffer[256];
1451: PetscViewer viewer;
1452: PetscViewerFormat format;
1453: PetscInt adjDim, adjCodim, overlap, minRefinement, initRefinement, maxRefinement, grade;
1454: PetscReal weightsFactor;
1455: DMForestAdaptivityStrategy adaptStrategy;
1456: PetscErrorCode ierr;
1460: DMForestGetTopology(dm, &oldTopo);
1461: PetscOptionsHead(PetscOptionsObject,"DMForest Options");
1462: PetscOptionsString("-dm_forest_topology","the topology of the forest's base mesh","DMForestSetTopology",oldTopo,stringBuffer,sizeof(stringBuffer),&flg1);
1463: PetscOptionsViewer("-dm_forest_base_dm","load the base DM from a viewer specification","DMForestSetBaseDM",&viewer,&format,&flg2);
1464: PetscOptionsViewer("-dm_forest_coarse_forest","load the coarse forest from a viewer specification","DMForestSetCoarseForest",&viewer,&format,&flg3);
1465: PetscOptionsViewer("-dm_forest_fine_forest","load the fine forest from a viewer specification","DMForestSetFineForest",&viewer,&format,&flg4);
1466: 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}");
1467: if (flg1) {
1468: DMForestSetTopology(dm,(DMForestTopology)stringBuffer);
1469: DMForestSetBaseDM(dm,NULL);
1470: DMForestSetAdaptivityForest(dm,NULL);
1471: }
1472: if (flg2) {
1473: DM base;
1475: DMCreate(PetscObjectComm((PetscObject)dm),&base);
1476: PetscViewerPushFormat(viewer,format);
1477: DMLoad(base,viewer);
1478: PetscViewerDestroy(&viewer);
1479: DMForestSetBaseDM(dm,base);
1480: DMDestroy(&base);
1481: DMForestSetTopology(dm,NULL);
1482: DMForestSetAdaptivityForest(dm,NULL);
1483: }
1484: if (flg3) {
1485: DM coarse;
1487: DMCreate(PetscObjectComm((PetscObject)dm),&coarse);
1488: PetscViewerPushFormat(viewer,format);
1489: DMLoad(coarse,viewer);
1490: PetscViewerDestroy(&viewer);
1491: DMForestSetAdaptivityForest(dm,coarse);
1492: DMDestroy(&coarse);
1493: DMForestSetTopology(dm,NULL);
1494: DMForestSetBaseDM(dm,NULL);
1495: }
1496: if (flg4) {
1497: DM fine;
1499: DMCreate(PetscObjectComm((PetscObject)dm),&fine);
1500: PetscViewerPushFormat(viewer,format);
1501: DMLoad(fine,viewer);
1502: PetscViewerDestroy(&viewer);
1503: DMForestSetAdaptivityForest(dm,fine);
1504: DMDestroy(&fine);
1505: DMForestSetTopology(dm,NULL);
1506: DMForestSetBaseDM(dm,NULL);
1507: }
1508: DMForestGetAdjacencyDimension(dm,&adjDim);
1509: PetscOptionsBoundedInt("-dm_forest_adjacency_dimension","set the dimension of points that define adjacency in the forest","DMForestSetAdjacencyDimension",adjDim,&adjDim,&flg,0);
1510: if (flg) {
1511: DMForestSetAdjacencyDimension(dm,adjDim);
1512: } else {
1513: DMForestGetAdjacencyCodimension(dm,&adjCodim);
1514: PetscOptionsBoundedInt("-dm_forest_adjacency_codimension","set the codimension of points that define adjacency in the forest","DMForestSetAdjacencyCodimension",adjCodim,&adjCodim,&flg,1);
1515: if (flg) {
1516: DMForestSetAdjacencyCodimension(dm,adjCodim);
1517: }
1518: }
1519: DMForestGetPartitionOverlap(dm,&overlap);
1520: PetscOptionsBoundedInt("-dm_forest_partition_overlap","set the degree of partition overlap","DMForestSetPartitionOverlap",overlap,&overlap,&flg,0);
1521: if (flg) {
1522: DMForestSetPartitionOverlap(dm,overlap);
1523: }
1524: #if 0
1525: 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);
1526: if (flg) {
1527: DMForestSetMinimumRefinement(dm,minRefinement);
1528: DMForestSetInitialRefinement(dm,minRefinement);
1529: }
1530: PetscOptionsBoundedInt("-dm_refine_hierarchy","equivalent to -dm_forest_set_minimum_refinement 0 and -dm_forest_set_initial_refinement",NULL,initRefinement,&initRefinement,&flg,0);
1531: if (flg) {
1532: DMForestSetMinimumRefinement(dm,0);
1533: DMForestSetInitialRefinement(dm,initRefinement);
1534: }
1535: #endif
1536: DMForestGetMinimumRefinement(dm,&minRefinement);
1537: PetscOptionsBoundedInt("-dm_forest_minimum_refinement","set the minimum level of refinement in the forest","DMForestSetMinimumRefinement",minRefinement,&minRefinement,&flg,0);
1538: if (flg) {
1539: DMForestSetMinimumRefinement(dm,minRefinement);
1540: }
1541: DMForestGetInitialRefinement(dm,&initRefinement);
1542: PetscOptionsBoundedInt("-dm_forest_initial_refinement","set the initial level of refinement in the forest","DMForestSetInitialRefinement",initRefinement,&initRefinement,&flg,0);
1543: if (flg) {
1544: DMForestSetInitialRefinement(dm,initRefinement);
1545: }
1546: DMForestGetMaximumRefinement(dm,&maxRefinement);
1547: PetscOptionsBoundedInt("-dm_forest_maximum_refinement","set the maximum level of refinement in the forest","DMForestSetMaximumRefinement",maxRefinement,&maxRefinement,&flg,0);
1548: if (flg) {
1549: DMForestSetMaximumRefinement(dm,maxRefinement);
1550: }
1551: DMForestGetAdaptivityStrategy(dm,&adaptStrategy);
1552: PetscOptionsString("-dm_forest_adaptivity_strategy","the forest's adaptivity-flag resolution strategy","DMForestSetAdaptivityStrategy",adaptStrategy,stringBuffer,sizeof(stringBuffer),&flg);
1553: if (flg) {
1554: DMForestSetAdaptivityStrategy(dm,(DMForestAdaptivityStrategy)stringBuffer);
1555: }
1556: DMForestGetGradeFactor(dm,&grade);
1557: PetscOptionsBoundedInt("-dm_forest_grade_factor","grade factor between neighboring cells","DMForestSetGradeFactor",grade,&grade,&flg,0);
1558: if (flg) {
1559: DMForestSetGradeFactor(dm,grade);
1560: }
1561: DMForestGetCellWeightFactor(dm,&weightsFactor);
1562: PetscOptionsReal("-dm_forest_cell_weight_factor","multiplying weight factor for cell refinement","DMForestSetCellWeightFactor",weightsFactor,&weightsFactor,&flg);
1563: if (flg) {
1564: DMForestSetCellWeightFactor(dm,weightsFactor);
1565: }
1566: PetscOptionsTail();
1567: return(0);
1568: }
1570: PetscErrorCode DMCreateSubDM_Forest(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm)
1571: {
1575: if (subdm) {DMClone(dm, subdm);}
1576: DMCreateSectionSubDM(dm, numFields, fields, is, subdm);
1577: return(0);
1578: }
1580: PetscErrorCode DMRefine_Forest(DM dm, MPI_Comm comm, DM *dmRefined)
1581: {
1582: DMLabel refine;
1583: DM fineDM;
1587: DMGetFineDM(dm,&fineDM);
1588: if (fineDM) {
1589: PetscObjectReference((PetscObject)fineDM);
1590: *dmRefined = fineDM;
1591: return(0);
1592: }
1593: DMForestTemplate(dm,comm,dmRefined);
1594: DMGetLabel(dm,"refine",&refine);
1595: if (!refine) {
1596: DMLabelCreate(PETSC_COMM_SELF, "refine",&refine);
1597: DMLabelSetDefaultValue(refine,DM_ADAPT_REFINE);
1598: } else {
1599: PetscObjectReference((PetscObject) refine);
1600: }
1601: DMForestSetAdaptivityLabel(*dmRefined,refine);
1602: DMLabelDestroy(&refine);
1603: return(0);
1604: }
1606: PetscErrorCode DMCoarsen_Forest(DM dm, MPI_Comm comm, DM *dmCoarsened)
1607: {
1608: DMLabel coarsen;
1609: DM coarseDM;
1613: {
1614: PetscMPIInt mpiComparison;
1615: MPI_Comm dmcomm = PetscObjectComm((PetscObject)dm);
1617: MPI_Comm_compare(comm, dmcomm, &mpiComparison);
1618: if (mpiComparison != MPI_IDENT && mpiComparison != MPI_CONGRUENT) SETERRQ(dmcomm,PETSC_ERR_SUP,"No support for different communicators yet");
1619: }
1620: DMGetCoarseDM(dm,&coarseDM);
1621: if (coarseDM) {
1622: PetscObjectReference((PetscObject)coarseDM);
1623: *dmCoarsened = coarseDM;
1624: return(0);
1625: }
1626: DMForestTemplate(dm,comm,dmCoarsened);
1627: DMForestSetAdaptivityPurpose(*dmCoarsened,DM_ADAPT_COARSEN);
1628: DMGetLabel(dm,"coarsen",&coarsen);
1629: if (!coarsen) {
1630: DMLabelCreate(PETSC_COMM_SELF, "coarsen",&coarsen);
1631: DMLabelSetDefaultValue(coarsen,DM_ADAPT_COARSEN);
1632: } else {
1633: PetscObjectReference((PetscObject) coarsen);
1634: }
1635: DMForestSetAdaptivityLabel(*dmCoarsened,coarsen);
1636: DMLabelDestroy(&coarsen);
1637: return(0);
1638: }
1640: static PetscErrorCode DMAdaptLabel_Forest(DM dm, DMLabel label, DM *adaptedDM)
1641: {
1642: PetscBool success;
1646: DMForestTemplate(dm,PetscObjectComm((PetscObject)dm),adaptedDM);
1647: DMForestSetAdaptivityLabel(*adaptedDM,label);
1648: DMSetUp(*adaptedDM);
1649: DMForestGetAdaptivitySuccess(*adaptedDM,&success);
1650: if (!success) {
1651: DMDestroy(adaptedDM);
1652: *adaptedDM = NULL;
1653: }
1654: return(0);
1655: }
1657: static PetscErrorCode DMInitialize_Forest(DM dm)
1658: {
1662: PetscMemzero(dm->ops,sizeof(*(dm->ops)));
1664: dm->ops->clone = DMClone_Forest;
1665: dm->ops->setfromoptions = DMSetFromOptions_Forest;
1666: dm->ops->destroy = DMDestroy_Forest;
1667: dm->ops->createsubdm = DMCreateSubDM_Forest;
1668: dm->ops->refine = DMRefine_Forest;
1669: dm->ops->coarsen = DMCoarsen_Forest;
1670: dm->ops->adaptlabel = DMAdaptLabel_Forest;
1671: return(0);
1672: }
1674: /*MC
1676: DMFOREST = "forest" - A DM object that encapsulates a hierarchically refined mesh. Forests usually have a base DM
1677: (see DMForestGetBaseDM()), from which it is refined. The refinement and partitioning of forests is considered
1678: immutable after DMSetUp() is called. To adapt a mesh, one should call DMForestTemplate() to create a new mesh that
1679: will default to being identical to it, specify how that mesh should differ, and then calling DMSetUp() on the new
1680: mesh.
1682: To specify that a mesh should be refined or coarsened from the previous mesh, a label should be defined on the
1683: previous mesh whose values indicate which cells should be refined (DM_ADAPT_REFINE) or coarsened (DM_ADAPT_COARSEN)
1684: and how (subtypes are free to allow additional values for things like anisotropic refinement). The label should be
1685: given to the *new* mesh with DMForestSetAdaptivityLabel().
1687: Level: advanced
1689: .seealso: DMType, DMCreate(), DMSetType(), DMForestGetBaseDM(), DMForestSetBaseDM(), DMForestTemplate(), DMForestSetAdaptivityLabel()
1690: M*/
1692: PETSC_EXTERN PetscErrorCode DMCreate_Forest(DM dm)
1693: {
1694: DM_Forest *forest;
1699: PetscNewLog(dm,&forest);
1700: dm->dim = 0;
1701: dm->data = forest;
1702: forest->refct = 1;
1703: forest->data = NULL;
1704: forest->topology = NULL;
1705: forest->adapt = NULL;
1706: forest->base = NULL;
1707: forest->adaptPurpose = DM_ADAPT_DETERMINE;
1708: forest->adjDim = PETSC_DEFAULT;
1709: forest->overlap = PETSC_DEFAULT;
1710: forest->minRefinement = PETSC_DEFAULT;
1711: forest->maxRefinement = PETSC_DEFAULT;
1712: forest->initRefinement = PETSC_DEFAULT;
1713: forest->cStart = PETSC_DETERMINE;
1714: forest->cEnd = PETSC_DETERMINE;
1715: forest->cellSF = NULL;
1716: forest->adaptLabel = NULL;
1717: forest->gradeFactor = 2;
1718: forest->cellWeights = NULL;
1719: forest->cellWeightsCopyMode = PETSC_USE_POINTER;
1720: forest->weightsFactor = 1.;
1721: forest->weightCapacity = 1.;
1722: DMForestSetAdaptivityStrategy(dm,DMFORESTADAPTALL);
1723: DMInitialize_Forest(dm);
1724: return(0);
1725: }