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