1 /* 2 DMForest, for parallel, hierarchically refined, distributed mesh problems. 3 */ 4 #if !defined(__PETSCDMFOREST_H) 5 #define __PETSCDMFOREST_H 6 7 #include <petscdm.h> 8 9 /*J 10 DMForestTopology - String with the name of a PETSc DMForest base mesh topology 11 12 Level: beginner 13 14 .seealso: DMForestSetType(), DMForest 15 J*/ 16 typedef const char* DMForestTopology; 17 18 PETSC_EXTERN PetscErrorCode DMForestSetTopology(DM, DMForestTopology); 19 PETSC_EXTERN PetscErrorCode DMForestGetTopology(DM, DMForestTopology *); 20 21 /* this is the coarsest possible forest: can be any DM which we can 22 * convert to a DMForest */ 23 PETSC_EXTERN PetscErrorCode DMForestSetBaseDM(DM, DM); 24 PETSC_EXTERN PetscErrorCode DMForestGetBaseDM(DM, DM *); 25 26 /* these are forests from which we can adapt */ 27 PETSC_EXTERN PetscErrorCode DMForestSetCoarseForest(DM, DM); 28 PETSC_EXTERN PetscErrorCode DMForestGetCoarseForest(DM, DM *); 29 PETSC_EXTERN PetscErrorCode DMForestSetFineForest(DM, DM); 30 PETSC_EXTERN PetscErrorCode DMForestGetFineForest(DM, DM *); 31 32 /* what we consider adjacent, for the purposes of cell grading, overlap, etc. */ 33 PETSC_EXTERN PetscErrorCode DMForestSetAdjacencyDimension(DM, PetscInt); 34 PETSC_EXTERN PetscErrorCode DMForestGetAdjacencyDimension(DM, PetscInt *); 35 PETSC_EXTERN PetscErrorCode DMForestSetAdjacencyCodimension(DM, PetscInt); 36 PETSC_EXTERN PetscErrorCode DMForestGetAdjacencyCodimension(DM, PetscInt *); 37 38 PETSC_EXTERN PetscErrorCode DMForestSetPartitionOverlap(DM, PetscInt); 39 PETSC_EXTERN PetscErrorCode DMForestGetPartitionOverlap(DM, PetscInt *); 40 41 PETSC_EXTERN PetscErrorCode DMForestSetMinimumRefinement(DM, PetscInt); 42 PETSC_EXTERN PetscErrorCode DMForestGetMinimumRefinement(DM, PetscInt *); 43 44 PETSC_EXTERN PetscErrorCode DMForestSetMaximumRefinement(DM, PetscInt); 45 PETSC_EXTERN PetscErrorCode DMForestGetMaximumRefinement(DM, PetscInt *); 46 47 PETSC_EXTERN PetscErrorCode DMForestSetInitialRefinement(DM, PetscInt); 48 PETSC_EXTERN PetscErrorCode DMForestGetInitialRefinement(DM, PetscInt *); 49 50 PETSC_EXTERN PetscErrorCode DMForestGetCellChart(DM, PetscInt *, PetscInt *); 51 PETSC_EXTERN PetscErrorCode DMForestGetCellSF(DM, PetscSF *); 52 53 /* reserve some adaptivity types */ 54 enum {DM_FOREST_KEEP = 0, 55 DM_FOREST_REFINE, 56 DM_FOREST_COARSEN, 57 DM_FOREST_RESERVED_ADAPTIVITY_COUNT}; 58 59 /* flag each cell with an adaptivity count: should match the cell section */ 60 PETSC_EXTERN PetscErrorCode DMForestSetCellAdaptivityMarkers(DM, PetscInt[], PetscCopyMode); 61 PETSC_EXTERN PetscErrorCode DMForestGetCellAdaptivityMarkers(DM, PetscInt *[]); 62 63 /*J 64 DMForestAdaptivityStrategy - String with the name of a PETSc DMForest adaptivity strategy 65 66 Level: intermediate 67 68 .seealso: DMForestSetType(), DMForest 69 J*/ 70 typedef const char* DMForestAdaptivityStrategy; 71 #define DMFORESTADAPTALL "all" 72 #define DMFORESTADAPTANY "any" 73 74 /* how to combine: -flags from multiple processes, 75 * -coarsen flags from multiple children 76 */ 77 PETSC_EXTERN PetscErrorCode DMForestSetAdaptivityStrategy(DM, DMForestAdaptivityStrategy); 78 PETSC_EXTERN PetscErrorCode DMForestGetAdaptivityStrategy(DM, DMForestAdaptivityStrategy *); 79 80 /* for a quadtree/octree mesh, this is the x:1 condition: 1 indicates a uniform mesh, 81 * 2 indicates typical 2:1, 82 */ 83 PETSC_EXTERN PetscErrorCode DMForestSetGradeFactor(DM, PetscInt); 84 PETSC_EXTERN PetscErrorCode DMForestGetGradeFactor(DM, PetscInt *); 85 86 /* weights for repartitioning */ 87 PETSC_EXTERN PetscErrorCode DMForestSetCellWeights(DM, PetscReal[], PetscCopyMode); 88 PETSC_EXTERN PetscErrorCode DMForestGetCellWeights(DM, PetscReal *[]); 89 90 /* weight multiplier for refinement level: useful for sub-cycling time steps */ 91 PETSC_EXTERN PetscErrorCode DMForestSetCellWeightFactor(DM, PetscReal); 92 PETSC_EXTERN PetscErrorCode DMForestGetCellWeightFactor(DM, PetscReal *); 93 94 /* this process's capacity when redistributing the cells */ 95 PETSC_EXTERN PetscErrorCode DMForestSetWeightCapacity(DM, PetscReal); 96 PETSC_EXTERN PetscErrorCode DMForestGetWeightCapacity(DM, PetscReal *); 97 98 PETSC_EXTERN PetscErrorCode DMForestSetFromOptions(DM); 99 PETSC_EXTERN PetscErrorCode DMForestSetUp(DM); 100 101 PETSC_EXTERN PetscErrorCode DMForestGetFineProjector(DM,Mat *); 102 PETSC_EXTERN PetscErrorCode DMForestGetCoarseRestrictor(DM,Mat *); 103 104 #endif 105