1 #ifndef PETSCDMTYPES_H 2 #define PETSCDMTYPES_H 3 4 /* SUBMANSEC = DM */ 5 6 /*S 7 DM - Abstract PETSc object that manages an abstract grid object and its interactions with the algebraic solvers 8 9 Level: intermediate 10 11 .seealso: `DMType`, `DMDGetType()`, `DMCompositeCreate()`, `DMDACreate()`, `DMSetType()`, `DMType`, `DMDA`, `DMPLEX` 12 S*/ 13 typedef struct _p_DM *DM; 14 15 /*E 16 DMBoundaryType - Describes the choice for fill of ghost cells on physical domain boundaries. 17 18 Level: beginner 19 20 A boundary may be of type `DM_BOUNDARY_NONE` (no ghost nodes), `DM_BOUNDARY_GHOSTED` (ghost vertices/cells 21 exist but aren't filled; you can put values into them and then apply a stencil that uses those ghost locations), 22 `DM_BOUNDARY_MIRROR` (the ghost value is the same as the value 1 grid point in; that is, the 0th grid point in the real mesh acts like a mirror to define the ghost point value; 23 not yet implemented for 3d), `DM_BOUNDARY_PERIODIC` (ghost vertices/cells filled by the opposite 24 edge of the domain), or `DM_BOUNDARY_TWIST` (like periodic, only glued backwards like a Mobius strip). 25 26 Notes: 27 This is information for the boundary of the __PHYSICAL__ domain. It has nothing to do with boundaries between 28 processes. That width is always determined by the stencil width; see `DMDASetStencilWidth()`. 29 30 If the physical grid points have values 0 1 2 3 with `DM_BOUNDARY_MIRROR` then the local vector with ghost points has the values 1 0 1 2 3 2 . 31 32 Developer Notes: 33 Should` DM_BOUNDARY_MIRROR` have the same meaning with DMDA_Q0, that is a staggered grid? In that case should the ghost point have the same value 34 as the 0th grid point where the physical boundary serves as the mirror? 35 36 References: 37 . * - https://scicomp.stackexchange.com/questions/5355/writing-the-poisson-equation-finite-difference-matrix-with-neumann-boundary-cond 38 39 .seealso: `DMDASetBoundaryType()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`, `DMDACreate()` 40 E*/ 41 typedef enum { 42 DM_BOUNDARY_NONE, 43 DM_BOUNDARY_GHOSTED, 44 DM_BOUNDARY_MIRROR, 45 DM_BOUNDARY_PERIODIC, 46 DM_BOUNDARY_TWIST 47 } DMBoundaryType; 48 /*E 49 DMBoundaryConditionType - indicates what type of boundary condition is to be imposed 50 51 Note: This flag indicates the type of function which will define the condition: 52 $ DM_BC_ESSENTIAL - A Dirichlet condition using a function of the coordinates 53 $ DM_BC_ESSENTIAL_FIELD - A Dirichlet condition using a function of the coordinates and auxiliary field data 54 $ DM_BC_ESSENTIAL_BD_FIELD - A Dirichlet condition using a function of the coordinates, facet normal, and auxiliary field data 55 $ DM_BC_NATURAL - A Neumann condition using a function of the coordinates 56 $ DM_BC_NATURAL_FIELD - A Neumann condition using a function of the coordinates and auxiliary field data 57 $ DM_BC_NATURAL_RIEMANN - A flux condition which determines the state in ghost cells 58 The user can check whether a boundary condition is essential using (type & DM_BC_ESSENTIAL), and similarly for 59 natural conditions (type & DM_BC_NATURAL) 60 61 Level: beginner 62 63 .seealso: `DMAddBoundary()`, `DSAddBoundary()`, `DSGetBoundary()` 64 E*/ 65 typedef enum { 66 DM_BC_ESSENTIAL = 1, 67 DM_BC_ESSENTIAL_FIELD = 5, 68 DM_BC_NATURAL = 2, 69 DM_BC_NATURAL_FIELD = 6, 70 DM_BC_ESSENTIAL_BD_FIELD = 9, 71 DM_BC_NATURAL_RIEMANN = 10 72 } DMBoundaryConditionType; 73 74 /*E 75 DMPointLocationType - Describes the method to handle point location failure 76 77 Level: beginner 78 79 If a search using `DM_POINTLOCATION_NONE` fails, the failure is signaled with a negative cell number. On the 80 other hand, if `DM_POINTLOCATION_NEAREST` is used, on failure, the (approximate) nearest point in the mesh is 81 used, replacing the given point in the input vector. `DM_POINTLOCATION_REMOVE` returns values only for points 82 which were located. 83 84 .seealso: `DMLocatePoints()` 85 E*/ 86 typedef enum { 87 DM_POINTLOCATION_NONE, 88 DM_POINTLOCATION_NEAREST, 89 DM_POINTLOCATION_REMOVE 90 } DMPointLocationType; 91 92 /*E 93 DMAdaptationStrategy - Describes the strategy used for adaptive solves 94 95 Level: beginner 96 97 DM_ADAPTATION_INITIAL will refine a mesh based on an initial guess. DM_ADAPTATION_SEQUENTIAL will refine the 98 mesh based on a sequence of solves, much like grid sequencing. DM_ADAPTATION_MULTILEVEL will use the sequence 99 of constructed meshes in a multilevel solve, much like the Systematic Upscaling of Brandt. 100 101 .seealso: `DMAdaptorSolve()` 102 E*/ 103 typedef enum { 104 DM_ADAPTATION_INITIAL, 105 DM_ADAPTATION_SEQUENTIAL, 106 DM_ADAPTATION_MULTILEVEL 107 } DMAdaptationStrategy; 108 109 /*E 110 DMAdaptationCriterion - Describes the test used to decide whether to coarsen or refine parts of the mesh 111 112 Level: beginner 113 114 `DM_ADAPTATION_REFINE` will uniformly refine a mesh, much like grid sequencing. `DM_ADAPTATION_LABEL` will adapt 115 the mesh based upon a label of the cells filled with `DMAdaptFlag` markers. `DM_ADAPTATION_METRIC` will try to 116 mesh the manifold described by the input metric tensor uniformly. PETSc can also construct such a metric based 117 upon an input primal or a gradient field. 118 119 .seealso: `DMAdaptorSolve()` 120 E*/ 121 typedef enum { 122 DM_ADAPTATION_NONE, 123 DM_ADAPTATION_REFINE, 124 DM_ADAPTATION_LABEL, 125 DM_ADAPTATION_METRIC 126 } DMAdaptationCriterion; 127 128 /*E 129 DMAdaptFlag - Marker in the label prescribing adaptation 130 131 Level: beginner 132 133 .seealso: `DMAdaptLabel()` 134 E*/ 135 typedef enum { 136 DM_ADAPT_DETERMINE = PETSC_DETERMINE, 137 DM_ADAPT_KEEP = 0, 138 DM_ADAPT_REFINE, 139 DM_ADAPT_COARSEN, 140 DM_ADAPT_COARSEN_LAST, 141 DM_ADAPT_RESERVED_COUNT 142 } DMAdaptFlag; 143 144 /*E 145 DMDirection - Indicates a coordinate direction 146 147 Level: beginner 148 149 .seealso: `DMDAGetRay()`, `DMDAGetProcessorSubset()`, `DMPlexShearGeometry()` 150 E*/ 151 typedef enum { 152 DM_X, 153 DM_Y, 154 DM_Z 155 } DMDirection; 156 157 /*E 158 DMEnclosureType - The type of enclosure relation between one `DM` and another 159 160 Level: beginner 161 162 Notes: 163 For example, one `DM` dmA may be the boundary of another dmB, in which case it would be labeled `DM_ENC_SUBMESH`. 164 165 If the situation is reversed, and dmA has boundary dmB, it would be labeled `DM_ENC_SUPERMESH`. 166 167 Likewise, if dmA was a subregion of dmB, it would be labeled `DM_ENC_SUBMESH`. 168 169 If no relation can be determined, `DM_ENC_NONE` is used. 170 171 If a relation is not yet known, then `DM_ENC_UNKNOWN` is used. 172 173 .seealso: `DMGetEnclosureRelation()` 174 E*/ 175 typedef enum { 176 DM_ENC_EQUALITY, 177 DM_ENC_SUPERMESH, 178 DM_ENC_SUBMESH, 179 DM_ENC_NONE, 180 DM_ENC_UNKNOWN 181 } DMEnclosureType; 182 183 /*E 184 DMPolytopeType - This describes the polytope represented by each cell. 185 186 Level: beginner 187 188 While most operations only need the topology information in the `DMPLEX`, we must sometimes have the 189 user specify a polytope. For instance, when interpolating from a cell-vertex mesh, the type of 190 polytope can be ambiguous. Also, `DMPLEX` allows different symmetries of a prism cell with the same 191 constituent points. Normally these types are autoamtically inferred and the user does not specify 192 them. 193 194 .seealso: `DMPlexComputeCellTypes()` 195 E*/ 196 typedef enum { 197 DM_POLYTOPE_POINT, 198 DM_POLYTOPE_SEGMENT, 199 DM_POLYTOPE_POINT_PRISM_TENSOR, 200 DM_POLYTOPE_TRIANGLE, 201 DM_POLYTOPE_QUADRILATERAL, 202 DM_POLYTOPE_SEG_PRISM_TENSOR, 203 DM_POLYTOPE_TETRAHEDRON, 204 DM_POLYTOPE_HEXAHEDRON, 205 DM_POLYTOPE_TRI_PRISM, 206 DM_POLYTOPE_TRI_PRISM_TENSOR, 207 DM_POLYTOPE_QUAD_PRISM_TENSOR, 208 DM_POLYTOPE_PYRAMID, 209 DM_POLYTOPE_FV_GHOST, 210 DM_POLYTOPE_INTERIOR_GHOST, 211 DM_POLYTOPE_UNKNOWN, 212 DM_NUM_POLYTOPES 213 } DMPolytopeType; 214 PETSC_EXTERN const char *const DMPolytopeTypes[]; 215 216 /*E 217 PetscUnit - The seven fundamental SI units 218 219 Level: beginner 220 221 .seealso: `DMPlexGetScale()`, `DMPlexSetScale()` 222 E*/ 223 typedef enum { 224 PETSC_UNIT_LENGTH, 225 PETSC_UNIT_MASS, 226 PETSC_UNIT_TIME, 227 PETSC_UNIT_CURRENT, 228 PETSC_UNIT_TEMPERATURE, 229 PETSC_UNIT_AMOUNT, 230 PETSC_UNIT_LUMINOSITY, 231 NUM_PETSC_UNITS 232 } PetscUnit; 233 234 /*S 235 DMField - PETSc object for defining a field on a mesh topology 236 237 Level: intermediate 238 S*/ 239 typedef struct _p_DMField *DMField; 240 241 /*S 242 DMUniversalLabel - A label that encodes a set of `DMLabel`s, bijectively 243 244 Level: developer 245 S*/ 246 typedef struct _p_UniversalLabel *DMUniversalLabel; 247 248 typedef struct _n_DMGeneratorFunctionList *DMGeneratorFunctionList; 249 250 #endif 251