1 #pragma once 2 3 /* SUBMANSEC = DM */ 4 5 /*S 6 DM - Abstract PETSc object that manages an abstract grid-like object and its interactions with the algebraic solvers 7 8 Level: intermediate 9 10 Note: 11 `DM` is an orphan initialism or orphan acronym, the letters have no meaning and never did. 12 13 .seealso: [](ch_dmbase), `DMType`, `DMGetType()`, `DMCompositeCreate()`, `DMDACreate()`, `DMSetType()`, `DMType`, `DMDA`, `DMPLEX` 14 S*/ 15 typedef struct _p_DM *DM; 16 17 /*E 18 DMBoundaryType - Describes the choice for the filling of ghost cells on physical domain boundaries. 19 20 Values: 21 + `DM_BOUNDARY_NONE` - no ghost nodes 22 . `DM_BOUNDARY_GHOSTED` - ghost vertices/cells exist but aren't filled; you can put values into them and then apply a stencil that uses those ghost locations 23 . `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 24 the ghost point value; not yet implemented for 3d 25 . `DM_BOUNDARY_PERIODIC` - ghost vertices/cells filled by the opposite edge of the domain 26 - `DM_BOUNDARY_TWIST` - like periodic, only glued backwards like a Mobius strip 27 28 Level: beginner 29 30 Notes: 31 This is information for the boundary of the __PHYSICAL__ domain. It has nothing to do with boundaries between 32 processes. That width is always determined by the stencil width; see `DMDASetStencilWidth()`. 33 34 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. 35 36 See <https://scicomp.stackexchange.com/questions/5355/writing-the-poisson-equation-finite-difference-matrix-with-neumann-boundary-cond> 37 38 Developer Note: 39 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 40 as the 0th grid point where the physical boundary serves as the mirror? 41 42 .seealso: [](ch_dmbase), `DM`, `DMDA`, `DMDASetBoundaryType()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`, `DMDACreate()` 43 E*/ 44 typedef enum { 45 DM_BOUNDARY_NONE, 46 DM_BOUNDARY_GHOSTED, 47 DM_BOUNDARY_MIRROR, 48 DM_BOUNDARY_PERIODIC, 49 DM_BOUNDARY_TWIST 50 } DMBoundaryType; 51 52 /*E 53 DMBoundaryConditionType - indicates what type of boundary condition is to be imposed 54 55 Values: 56 + `DM_BC_ESSENTIAL` - A Dirichlet condition using a function of the coordinates 57 . `DM_BC_ESSENTIAL_FIELD` - A Dirichlet condition using a function of the coordinates and auxiliary field data 58 . `DM_BC_ESSENTIAL_BD_FIELD` - A Dirichlet condition using a function of the coordinates, facet normal, and auxiliary field data 59 . `DM_BC_NATURAL` - A Neumann condition using a function of the coordinates 60 . `DM_BC_NATURAL_FIELD` - A Neumann condition using a function of the coordinates and auxiliary field data 61 . `DM_BC_NATURAL_RIEMANN` - A flux condition which determines the state in ghost cells 62 . `DM_BC_LOWER_BOUND` - A lower bound on the solution along a boundary 63 - `DM_BC_UPPER_BOUND` - An upper bound on the solution along a boundary 64 65 Level: beginner 66 67 Note: 68 The user can check whether a boundary condition is essential using (type & `DM_BC_ESSENTIAL`), and similarly for 69 natural conditions (type & `DM_BC_NATURAL`) 70 71 .seealso: [](ch_dmbase), `DM`, `DMAddBoundary()`, `DSAddBoundary()`, `DSGetBoundary()` 72 E*/ 73 typedef enum { 74 DM_BC_ESSENTIAL = 1, 75 DM_BC_ESSENTIAL_FIELD = 5, 76 DM_BC_NATURAL = 2, 77 DM_BC_NATURAL_FIELD = 6, 78 DM_BC_ESSENTIAL_BD_FIELD = 9, 79 DM_BC_NATURAL_RIEMANN = 10, 80 DM_BC_LOWER_BOUND = 4, 81 DM_BC_UPPER_BOUND = 8 82 } DMBoundaryConditionType; 83 84 /*E 85 DMPointLocationType - Describes the method to handle point location failure 86 87 Values: 88 + `DM_POINTLOCATION_NONE` - return a negative cell number 89 . `DM_POINTLOCATION_NEAREST` - the (approximate) nearest point in the mesh is used 90 - `DM_POINTLOCATION_REMOVE` - returns values only for points which were located 91 92 Level: intermediate 93 94 .seealso: [](ch_dmbase), `DM`, `DMLocatePoints()` 95 E*/ 96 typedef enum { 97 DM_POINTLOCATION_NONE, 98 DM_POINTLOCATION_NEAREST, 99 DM_POINTLOCATION_REMOVE 100 } DMPointLocationType; 101 102 /*E 103 DMBlockingType - Describes how to choose variable block sizes 104 105 Values: 106 + `DM_BLOCKING_TOPOLOGICAL_POINT` - select all fields at a topological point (cell center, at a face, etc) 107 - `DM_BLOCKING_FIELD_NODE` - using a separate block for each field at a topological point 108 109 Level: intermediate 110 111 Note: 112 When using `PCVPBJACOBI`, one can choose to block by topological point (all fields at a cell center, at a face, etc.) 113 or by field nodes (using number of components per field to identify "nodes"). Field nodes lead to smaller blocks, but 114 may converge more slowly. For example, a cubic Lagrange hexahedron will have one node at vertices, two at edges, four 115 at faces, and eight at cell centers. If using point blocking, the `PCVPBJACOBI` preconditioner will work with block 116 sizes up to 8 Lagrange nodes. For 5-component CFD, this produces matrices up to 40x40, which increases memory 117 footprint and may harm performance. With field node blocking, the maximum block size will correspond to one Lagrange node, 118 or 5x5 blocks for the CFD example. 119 120 .seealso: [](ch_dmbase), `PCVPBJACOBI`, `MatSetVariableBlockSizes()`, `DMSetBlockingType()` 121 E*/ 122 typedef enum { 123 DM_BLOCKING_TOPOLOGICAL_POINT, 124 DM_BLOCKING_FIELD_NODE 125 } DMBlockingType; 126 127 /*E 128 DMAdaptationStrategy - Describes the strategy used for adaptive solves 129 130 Values: 131 + `DM_ADAPTATION_INITIAL` - refine a mesh based on an initial guess 132 . `DM_ADAPTATION_SEQUENTIAL` - refine the mesh based on a sequence of solves, much like grid sequencing 133 - `DM_ADAPTATION_MULTILEVEL` - use the sequence of constructed meshes in a multilevel solve, much like the Systematic Upscaling of Brandt 134 135 Level: beginner 136 137 .seealso: [](ch_dmbase), `DM`, `DMAdaptor`, `DMAdaptationCriterion`, `DMAdaptorSolve()` 138 E*/ 139 typedef enum { 140 DM_ADAPTATION_INITIAL, 141 DM_ADAPTATION_SEQUENTIAL, 142 DM_ADAPTATION_MULTILEVEL 143 } DMAdaptationStrategy; 144 145 /*E 146 DMAdaptationCriterion - Describes the test used to decide whether to coarsen or refine parts of the mesh 147 148 Values: 149 + `DM_ADAPTATION_REFINE` - uniformly refine a mesh, much like grid sequencing 150 . `DM_ADAPTATION_LABEL` - adapt the mesh based upon a label of the cells filled with `DMAdaptFlag` markers. 151 . `DM_ADAPTATION_METRIC` - try to mesh the manifold described by the input metric tensor uniformly. PETSc can also construct such a metric based 152 upon an input primal or a gradient field. 153 - `DM_ADAPTATION_NONE` - do no adaptation 154 155 Level: beginner 156 157 .seealso: [](ch_dmbase), `DM`, `DMAdaptor`, `DMAdaptationStrategy`, `DMAdaptorSolve()` 158 E*/ 159 typedef enum { 160 DM_ADAPTATION_NONE, 161 DM_ADAPTATION_REFINE, 162 DM_ADAPTATION_LABEL, 163 DM_ADAPTATION_METRIC 164 } DMAdaptationCriterion; 165 PETSC_EXTERN const char *const DMAdaptationCriteria[]; 166 167 /*E 168 DMAdaptFlag - Marker in the label prescribing what adaptation to perform 169 170 Values: 171 + `DM_ADAPT_DETERMINE` - undocumented 172 . `DM_ADAPT_KEEP` - undocumented 173 . `DM_ADAPT_REFINE` - undocumented 174 . `DM_ADAPT_COARSEN` - undocumented 175 - `DM_ADAPT_COARSEN_LAST` - undocumented 176 177 Level: beginner 178 179 .seealso: [](ch_dmbase), `DM`, `DMAdaptor`, `DMAdaptationStrategy`, `DMAdaptationCriterion`, `DMAdaptorSolve()`, `DMAdaptLabel()` 180 E*/ 181 typedef enum { 182 DM_ADAPT_DETERMINE = PETSC_DETERMINE, 183 DM_ADAPT_KEEP = 0, 184 DM_ADAPT_REFINE = 1, 185 DM_ADAPT_COARSEN = 2, 186 DM_ADAPT_COARSEN_LAST = 3, 187 DM_ADAPT_RESERVED_COUNT = 4 188 } DMAdaptFlag; 189 190 /*E 191 DMDirection - Indicates a coordinate direction 192 193 Values: 194 + `DM_X` - the x coordinate direction 195 . `DM_Y` - the y coordinate direction 196 - `DM_Z` - the z coordinate direction 197 198 Level: beginner 199 200 .seealso: [](ch_dmbase), `DM`, `DMDA`, `DMDAGetRay()`, `DMDAGetProcessorSubset()`, `DMPlexShearGeometry()` 201 E*/ 202 typedef enum { 203 DM_X, 204 DM_Y, 205 DM_Z 206 } DMDirection; 207 208 /*E 209 DMEnclosureType - The type of enclosure relation between one `DM` and another 210 211 Values: 212 + `DM_ENC_SUBMESH` - the `DM` is the boundary of another `DM` 213 . `DM_ENC_SUPERMESH` - the `DM` has the boundary of another `DM` (the reverse situation to `DM_ENC_SUBMESH`) 214 . `DM_ENC_EQUALITY` - it is unknown what this means 215 . `DM_ENC_NONE` - no relationship can be determined 216 - `DM_ENC_UNKNOWN` - the relationship is unknown 217 218 Level: beginner 219 220 .seealso: [](ch_dmbase), `DM`, `DMGetEnclosureRelation()` 221 E*/ 222 typedef enum { 223 DM_ENC_EQUALITY, 224 DM_ENC_SUPERMESH, 225 DM_ENC_SUBMESH, 226 DM_ENC_NONE, 227 DM_ENC_UNKNOWN 228 } DMEnclosureType; 229 230 /*E 231 DMPolytopeType - This describes the polytope represented by each cell. 232 233 Level: beginner 234 235 While most operations only need the topology information in the `DMPLEX`, we must sometimes have the 236 user specify a polytope. For instance, when interpolating from a cell-vertex mesh, the type of 237 polytope can be ambiguous. Also, `DMPLEX` allows different symmetries of a prism cell with the same 238 constituent points. Normally these types are automatically inferred and the user does not specify 239 them. 240 241 .seealso: [](ch_dmbase), `DM`, `DMPlexComputeCellTypes()` 242 E*/ 243 typedef enum { 244 DM_POLYTOPE_POINT, 245 DM_POLYTOPE_SEGMENT, 246 DM_POLYTOPE_POINT_PRISM_TENSOR, 247 DM_POLYTOPE_TRIANGLE, 248 DM_POLYTOPE_QUADRILATERAL, 249 DM_POLYTOPE_SEG_PRISM_TENSOR, 250 DM_POLYTOPE_TETRAHEDRON, 251 DM_POLYTOPE_HEXAHEDRON, 252 DM_POLYTOPE_TRI_PRISM, 253 DM_POLYTOPE_TRI_PRISM_TENSOR, 254 DM_POLYTOPE_QUAD_PRISM_TENSOR, 255 DM_POLYTOPE_PYRAMID, 256 DM_POLYTOPE_FV_GHOST, 257 DM_POLYTOPE_INTERIOR_GHOST, 258 DM_POLYTOPE_UNKNOWN, 259 DM_POLYTOPE_UNKNOWN_CELL, 260 DM_POLYTOPE_UNKNOWN_FACE, 261 DM_NUM_POLYTOPES 262 } DMPolytopeType; 263 PETSC_EXTERN const char *const DMPolytopeTypes[]; 264 265 /*E 266 PetscUnit - The seven fundamental SI units 267 268 Level: beginner 269 270 .seealso: `DMPlexGetScale()`, `DMPlexSetScale()` 271 E*/ 272 typedef enum { 273 PETSC_UNIT_LENGTH, 274 PETSC_UNIT_MASS, 275 PETSC_UNIT_TIME, 276 PETSC_UNIT_CURRENT, 277 PETSC_UNIT_TEMPERATURE, 278 PETSC_UNIT_AMOUNT, 279 PETSC_UNIT_LUMINOSITY, 280 NUM_PETSC_UNITS 281 } PetscUnit; 282 283 /*E 284 DMReorderDefaultFlag - Flag indicating whether the `DM` should be reordered by default 285 286 Values: 287 + `DM_REORDER_DEFAULT_NOTSET` - Flag not set. 288 . `DM_REORDER_DEFAULT_FALSE` - Do not reorder by default. 289 - `DM_REORDER_DEFAULT_TRUE` - Reorder by default. 290 291 Level: intermediate 292 293 Developer Note: 294 Could be replaced with `PETSC_BOOL3` 295 296 .seealso: `DMPlexReorderSetDefault()`, `DMPlexReorderGetDefault()`, `DMPlexGetOrdering()`, `DMPlexPermute()` 297 E*/ 298 typedef enum { 299 DM_REORDER_DEFAULT_NOTSET = -1, 300 DM_REORDER_DEFAULT_FALSE = 0, 301 DM_REORDER_DEFAULT_TRUE = 1 302 } DMReorderDefaultFlag; 303 304 /*S 305 DMField - PETSc object for defining a field on a mesh topology 306 307 Level: intermediate 308 309 .seealso: [](ch_dmbase), `DM`, `DMUniversalLabel`, `DMLabelCreate()` 310 S*/ 311 typedef struct _p_DMField *DMField; 312 313 /*S 314 DMUniversalLabel - A label that encodes a set of `DMLabel`s, bijectively 315 316 Level: developer 317 318 .seealso: [](ch_dmbase), `DM`, `DMLabel`, `DMUniversalLabelCreate()` 319 S*/ 320 typedef struct _p_UniversalLabel *DMUniversalLabel; 321 322 typedef struct _PETSc_DMCEED *DMCeed; 323 324 typedef struct _n_DMGeneratorFunctionList *DMGeneratorFunctionList; 325