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 /* Just a name for the shape of the domain */ 19 PETSC_EXTERN PetscErrorCode DMForestSetTopology(DM, DMForestTopology); 20 PETSC_EXTERN PetscErrorCode DMForestGetTopology(DM, DMForestTopology *); 21 22 /* this is the coarsest possible forest: can be any DM which we can 23 * convert to a DMForest (right now: plex) */ 24 PETSC_EXTERN PetscErrorCode DMForestSetBaseDM(DM, DM); 25 PETSC_EXTERN PetscErrorCode DMForestGetBaseDM(DM, DM *); 26 PETSC_EXTERN PetscErrorCode DMForestSetBaseCoordinateMapping(DM, PetscErrorCode(*)(DM,PetscInt,PetscInt,const PetscReal[],PetscReal[],void*),void*); 27 PETSC_EXTERN PetscErrorCode DMForestGetBaseCoordinateMapping(DM, PetscErrorCode(**)(DM,PetscInt,PetscInt,const PetscReal[],PetscReal[],void*),void*); 28 29 /* this is the forest from which we adapt */ 30 PETSC_EXTERN PetscErrorCode DMForestSetAdaptivityForest(DM, DM); 31 PETSC_EXTERN PetscErrorCode DMForestGetAdaptivityForest(DM, DM *); 32 33 /* reserve some adaptivity types */ 34 enum {DM_FOREST_KEEP = 0, 35 DM_FOREST_REFINE, 36 DM_FOREST_COARSEN, 37 DM_FOREST_RESERVED_ADAPTIVITY_COUNT}; 38 39 typedef PetscInt DMForestAdaptivityPurpose; 40 PETSC_EXTERN PetscErrorCode DMForestSetAdaptivityPurpose(DM, DMForestAdaptivityPurpose); 41 PETSC_EXTERN PetscErrorCode DMForestGetAdaptivityPurpose(DM, DMForestAdaptivityPurpose*); 42 43 /* what we consider adjacent, for the purposes of cell grading, overlap, etc. */ 44 PETSC_EXTERN PetscErrorCode DMForestSetAdjacencyDimension(DM, PetscInt); 45 PETSC_EXTERN PetscErrorCode DMForestGetAdjacencyDimension(DM, PetscInt *); 46 PETSC_EXTERN PetscErrorCode DMForestSetAdjacencyCodimension(DM, PetscInt); 47 PETSC_EXTERN PetscErrorCode DMForestGetAdjacencyCodimension(DM, PetscInt *); 48 49 PETSC_EXTERN PetscErrorCode DMForestSetPartitionOverlap(DM, PetscInt); 50 PETSC_EXTERN PetscErrorCode DMForestGetPartitionOverlap(DM, PetscInt *); 51 52 PETSC_EXTERN PetscErrorCode DMForestSetMinimumRefinement(DM, PetscInt); 53 PETSC_EXTERN PetscErrorCode DMForestGetMinimumRefinement(DM, PetscInt *); 54 55 PETSC_EXTERN PetscErrorCode DMForestSetMaximumRefinement(DM, PetscInt); 56 PETSC_EXTERN PetscErrorCode DMForestGetMaximumRefinement(DM, PetscInt *); 57 58 PETSC_EXTERN PetscErrorCode DMForestSetInitialRefinement(DM, PetscInt); 59 PETSC_EXTERN PetscErrorCode DMForestGetInitialRefinement(DM, PetscInt *); 60 61 PETSC_EXTERN PetscErrorCode DMForestGetCellChart(DM, PetscInt *, PetscInt *); 62 PETSC_EXTERN PetscErrorCode DMForestGetCellSF(DM, PetscSF *); 63 64 65 /* flag each cell with an adaptivity count: should match the cell section */ 66 PETSC_EXTERN PetscErrorCode DMForestSetAdaptivityLabel(DM, const char *); 67 PETSC_EXTERN PetscErrorCode DMForestGetAdaptivityLabel(DM, const char **); 68 69 /*J 70 DMForestAdaptivityStrategy - String with the name of a PETSc DMForest adaptivity strategy 71 72 Level: intermediate 73 74 .seealso: DMForestSetType(), DMForest 75 J*/ 76 typedef const char* DMForestAdaptivityStrategy; 77 #define DMFORESTADAPTALL "all" 78 #define DMFORESTADAPTANY "any" 79 80 /* how to combine: -flags from multiple processes, 81 * -coarsen flags from multiple children 82 */ 83 PETSC_EXTERN PetscErrorCode DMForestSetAdaptivityStrategy(DM, DMForestAdaptivityStrategy); 84 PETSC_EXTERN PetscErrorCode DMForestGetAdaptivityStrategy(DM, DMForestAdaptivityStrategy *); 85 86 PETSC_EXTERN PetscErrorCode DMForestSetComputeAdaptivitySF(DM, PetscBool); 87 PETSC_EXTERN PetscErrorCode DMForestGetComputeAdaptivitySF(DM, PetscBool *); 88 89 PETSC_EXTERN PetscErrorCode DMForestGetAdaptivitySF(DM, PetscSF *, PetscSF *); 90 91 PETSC_EXTERN PetscErrorCode DMForestGetAdaptivitySuccess(DM, PetscBool *); 92 93 PETSC_EXTERN PetscErrorCode DMForestTransferVec(DM, Vec, DM, Vec, PetscBool, PetscReal); 94 95 /* for a quadtree/octree mesh, this is the x:1 condition: 1 indicates a uniform mesh, 96 * 2 indicates typical 2:1, 97 */ 98 PETSC_EXTERN PetscErrorCode DMForestSetGradeFactor(DM, PetscInt); 99 PETSC_EXTERN PetscErrorCode DMForestGetGradeFactor(DM, PetscInt *); 100 101 /* weights for repartitioning */ 102 PETSC_EXTERN PetscErrorCode DMForestSetCellWeights(DM, PetscReal[], PetscCopyMode); 103 PETSC_EXTERN PetscErrorCode DMForestGetCellWeights(DM, PetscReal *[]); 104 105 /* weight multiplier for refinement level: useful for sub-cycling time steps */ 106 PETSC_EXTERN PetscErrorCode DMForestSetCellWeightFactor(DM, PetscReal); 107 PETSC_EXTERN PetscErrorCode DMForestGetCellWeightFactor(DM, PetscReal *); 108 109 /* this process's capacity when redistributing the cells */ 110 PETSC_EXTERN PetscErrorCode DMForestSetWeightCapacity(DM, PetscReal); 111 PETSC_EXTERN PetscErrorCode DMForestGetWeightCapacity(DM, PetscReal *); 112 113 PETSC_EXTERN PetscErrorCode DMForestSetFromOptions(DM); 114 PETSC_EXTERN PetscErrorCode DMForestSetUp(DM); 115 116 /* convenience wrapper to common adaptivity pattern */ 117 PETSC_EXTERN PetscErrorCode DMForestAdaptLabel(DM, const char [], DM *); 118 119 PETSC_EXTERN PetscErrorCode DMForestGetFineProjector(DM,Mat *); 120 PETSC_EXTERN PetscErrorCode DMForestGetCoarseRestrictor(DM,Mat *); 121 122 /* miscellaneous */ 123 PETSC_EXTERN PetscErrorCode DMForestTemplate(DM,MPI_Comm,DM*); 124 125 /* type management */ 126 PETSC_EXTERN PetscErrorCode DMForestRegisterType(DMType); 127 PETSC_EXTERN PetscErrorCode DMIsForest(DM,PetscBool*); 128 129 /* p4est */ 130 PETSC_EXTERN PetscErrorCode DMP4estGetPartitionForCoarsening(DM,PetscBool *); 131 PETSC_EXTERN PetscErrorCode DMP4estSetPartitionForCoarsening(DM,PetscBool); 132 PETSC_EXTERN PetscErrorCode DMP8estGetPartitionForCoarsening(DM,PetscBool *); 133 PETSC_EXTERN PetscErrorCode DMP8estSetPartitionForCoarsening(DM,PetscBool); 134 135 #endif 136