1 #pragma once 2 3 #include <petscdmlabel.h> 4 5 /* MANSEC = DM */ 6 /* SUBMANSEC = DT */ 7 8 /*S 9 PetscDS - PETSc object that manages a discrete system, which is a set of discretizations + continuum equations from a `PetscWeakForm` 10 11 Level: intermediate 12 13 .seealso: `PetscDSCreate()`, `PetscDSSetType()`, `PetscDSType`, `PetscWeakForm`, `PetscFECreate()`, `PetscFVCreate()` 14 S*/ 15 typedef struct _p_PetscDS *PetscDS; 16 17 /*S 18 PetscWeakForm - PETSc object that manages a sets of pointwise functions defining a system of equations 19 20 Level: intermediate 21 22 .seealso: `PetscWeakFormCreate()`, `PetscDS`, `PetscFECreate()`, `PetscFVCreate()` 23 S*/ 24 typedef struct _p_PetscWeakForm *PetscWeakForm; 25 26 /*S 27 PetscFormKey - This key indicates how to use a set of pointwise functions defining part of a system of equations 28 29 The subdomain on which to integrate is specified by (label, value), the test function field by (field), and the 30 piece of the equation by (part). For example, LHS = 0 and RHS = 1 in IMEX methods. More pieces can be present for 31 operator splitting methods. 32 33 Level: intermediate 34 35 Note: 36 This is a struct, not a `PetscObject` 37 38 .seealso: `DMPlexSNESComputeResidualFEM()`, `DMPlexSNESComputeJacobianFEM()`, `DMPlexSNESComputeBoundaryFEM()` 39 S*/ 40 typedef struct { 41 DMLabel label; /* The (label, value) select a subdomain */ 42 PetscInt value; 43 PetscInt field; /* Selects the field for the test function */ 44 PetscInt part; /* Selects the equation part. For example, LHS = 0 and RHS = 1 in IMEX methods. More pieces can be present for operator splitting methods. */ 45 } PetscFormKey; 46 47 /*E 48 PetscWeakFormKind - The kind of weak form. The specific forms are given in the documentation for the integraton functions. 49 50 Values: 51 + OBJECTIVE - Objective form 52 . F0, F1 - Residual forms 53 . G0, G1, G2, G3 - Jacobian forms 54 . GP0, GP1, GP2, GP3 - Jacobian forms used to construct the preconditioner 55 . GT0, GT1, GT2, GT3 - Dynamic Jacobian matrix forms 56 . BDF0, BDF1 - Boundary Residual forms 57 . BDG0, BDG1, BDG2, BDG3 - Jacobian forms 58 . BDGP0, BDGP1, BDGP2, BDGP3 - Jacobian forms used to construct the preconditioner 59 . R - Riemann solver 60 - CEED - libCEED QFunction 61 62 Level: beginner 63 64 .seealso: `PetscWeakForm`, `PetscFEIntegrateResidual()`, `PetscFEIntegrateJacobian()`, `PetscFEIntegrateBdResidual()`, `PetscFEIntegrateBdJacobian()`, 65 `PetscFVIntegrateRHSFunction()`, `PetscWeakFormSetIndexResidual()`, `PetscWeakFormClearIndex()` 66 E*/ 67 typedef enum { 68 PETSC_WF_OBJECTIVE, 69 PETSC_WF_F0, 70 PETSC_WF_F1, 71 PETSC_WF_G0, 72 PETSC_WF_G1, 73 PETSC_WF_G2, 74 PETSC_WF_G3, 75 PETSC_WF_GP0, 76 PETSC_WF_GP1, 77 PETSC_WF_GP2, 78 PETSC_WF_GP3, 79 PETSC_WF_GT0, 80 PETSC_WF_GT1, 81 PETSC_WF_GT2, 82 PETSC_WF_GT3, 83 PETSC_WF_BDF0, 84 PETSC_WF_BDF1, 85 PETSC_WF_BDG0, 86 PETSC_WF_BDG1, 87 PETSC_WF_BDG2, 88 PETSC_WF_BDG3, 89 PETSC_WF_BDGP0, 90 PETSC_WF_BDGP1, 91 PETSC_WF_BDGP2, 92 PETSC_WF_BDGP3, 93 PETSC_WF_R, 94 PETSC_WF_CEED, 95 PETSC_NUM_WF 96 } PetscWeakFormKind; 97 PETSC_EXTERN const char *const PetscWeakFormKinds[]; 98 99 /*S 100 PetscPointFn - A prototype of a pointwise function that can be passed to, for example, `PetscDSSetObjective()` 101 102 Calling Sequence: 103 + dim - the coordinate dimension 104 . Nf - the number of fields 105 . NfAux - the number of auxiliary fields 106 . uOff - the offset into `u`[] and `u_t`[] for each field 107 . uOff_x - the offset into `u_x`[] for each field 108 . u - each field evaluated at the current point 109 . u_t - the time derivative of each field evaluated at the current point 110 . u_x - the gradient of each field evaluated at the current point 111 . aOff - the offset into `a`[] and `a_t`[] for each auxiliary field 112 . aOff_x - the offset into `a_x`[] for each auxiliary field 113 . a - each auxiliary field evaluated at the current point 114 . a_t - the time derivative of each auxiliary field evaluated at the current point 115 . a_x - the gradient of auxiliary each field evaluated at the current point 116 . t - current time 117 . x - coordinates of the current point 118 . numConstants - number of constant parameters 119 . constants - constant parameters 120 - obj - output values at the current point 121 122 Level: beginner 123 124 .seealso: `PetscPointFn`, `PetscDSSetObjective()`, `PetscDSGetObjective()`, PetscDSGetResidual()`, `PetscDSSetResidual()`, 125 `PetscDSGetRHSResidual()`, `PetscDSGetRHSResidual()`, `PetscDSSetUpdate()`, `PetscDSGetUpdate()`, `DMPlexSetCoordinateMap()` 126 S*/ 127 typedef void PetscPointFn(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal X[], PetscInt numConstants, const PetscScalar constants[], PetscScalar result[]); 128 129 /*S 130 PetscPointJacFn - A prototype of a pointwise function that can be passed to, for example, `PetscDSSetJacobian()` for computing Jacobians 131 132 Calling Sequence: 133 + dim - the coordinate dimension 134 . Nf - the number of fields 135 . NfAux - the number of auxiliary fields 136 . uOff - the offset into `u`[] and `u_t`[] for each field 137 . uOff_x - the offset into `u_x`[] for each field 138 . u - each field evaluated at the current point 139 . u_t - the time derivative of each field evaluated at the current point 140 . u_x - the gradient of each field evaluated at the current point 141 . aOff - the offset into `a`[] and `a_t`[] for each auxiliary field 142 . aOff_x - the offset into a_`x`[] for each auxiliary field 143 . a - each auxiliary field evaluated at the current point 144 . a_t - the time derivative of each auxiliary field evaluated at the current point 145 . a_x - the gradient of auxiliary each field evaluated at the current point 146 . t - current time 147 . u_tShift - the multiplier `a` for $dF/dU_t$ 148 . x - coordinates of the current point 149 . numConstants - number of constant parameters 150 . constants - constant parameters 151 - g - output values at the current point 152 153 Level: beginner 154 155 .seealso: `PetscPointFn`, `PetscDSSetJacobian()`, `PetscDSGetJacobian()`, PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobianPreconditioner()`, 156 `PetscDSSetDynamicJacobian()`, `PetscDSGetDynamicJacobian()` 157 S*/ 158 typedef void PetscPointJacFn(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g[]); 159 160 /*S 161 PetscBdPointFn - A prototype of a pointwise boundary function that can be passed to, for example, `PetscDSSetBdResidual()` 162 163 Calling Sequence: 164 + dim - the coordinate dimension 165 . Nf - the number of fields 166 . NfAux - the number of auxiliary fields 167 . uOff - the offset into `u`[] and `u_t`[] for each field 168 . uOff_x - the offset into `u_x`[] for each field 169 . u - each field evaluated at the current point 170 . u_t - the time derivative of each field evaluated at the current point 171 . u_x - the gradient of each field evaluated at the current point 172 . aOff - the offset into `a`[] and `a_t`[] for each auxiliary field 173 . aOff_x - the offset into `a_x`[] for each auxiliary field 174 . a - each auxiliary field evaluated at the current point 175 . a_t - the time derivative of each auxiliary field evaluated at the current point 176 . a_x - the gradient of auxiliary each field evaluated at the current point 177 . t - current time 178 . x - coordinates of the current point 179 . n - unit normal at the current point 180 . numConstants - number of constant parameters 181 . constants - constant parameters 182 - f - output values at the current point 183 184 Level: beginner 185 186 .seealso: `PetscPointFn`, `PetscDSSetBdResidual()`, `PetscDSGetBdResidual()`, `PetscDSSetObjective()`, `PetscDSGetObjective()`, PetscDSGetResidual()`, 187 `PetscDSGetRHSResidual()`, `PetscDSGetRHSResidual()`, `PetscDSSetUpdate()`, `PetscDSGetUpdate()`, `DMPlexSetCoordinateMap()`, 188 `PetscDSSetResidual()`, `PetscPointJacFn` 189 S*/ 190 typedef void PetscBdPointFn(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f[]); 191 192 /*S 193 PetscBdPointJacFn - A prototype of a pointwise boundary function that can be passed to, for example, `PetscDSSetBdJacobian()` 194 195 Calling Sequence: 196 + dim - the coordinate dimension 197 . Nf - the number of fields 198 . NfAux - the number of auxiliary fields 199 . uOff - the offset into `u`[] and `u_t`[] for each field 200 . uOff_x - the offset into `u_x`[] for each field 201 . u - each field evaluated at the current point 202 . u_t - the time derivative of each field evaluated at the current point 203 . u_x - the gradient of each field evaluated at the current point 204 . aOff - the offset into `a`[] and `a_t`[] for each auxiliary field 205 . aOff_x - the offset into `a_x`[] for each auxiliary field 206 . a - each auxiliary field evaluated at the current point 207 . a_t - the time derivative of each auxiliary field evaluated at the current point 208 . a_x - the gradient of auxiliary each field evaluated at the current point 209 . t - current time 210 . u_tShift - the multiplier `a` for $dF/dU_t$ 211 . x - coordinates of the current point 212 . n - normal at the current point 213 . numConstants - number of constant parameters 214 . constants - constant parameters 215 - g - output values at the current point 216 217 Level: beginner 218 219 .seealso: `PetscPointFn`, `PetscDSSetBdJacobian()`, PetscDSGetBdJacobian()`, `PetscDSSetBdJacobianPreconditioner()`, `PetscDSGetBdJacobianPreconditioner()`, 220 `PetscDSSetBdResidual()`, `PetscDSGetBdResidual()`, `PetscDSSetObjective()`, `PetscDSGetObjective()`, PetscDSGetResidual()`, 221 `PetscDSGetRHSResidual()`, `PetscDSGetRHSResidual()`, `PetscDSSetUpdate()`, `PetscDSGetUpdate()`, `DMPlexSetCoordinateMap()`, 222 `PetscDSSetResidual()`, `PetscPointJacFn` 223 S*/ 224 typedef void PetscBdPointJacFn(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]); 225 226 /*S 227 PetscPointExactSolutionFn - A prototype of a pointwise function that computes the exact solution to a PDE. Used with, for example, 228 `PetscDSSetExactSolution()` 229 230 Calling Sequence: 231 + dim - the coordinate dimension 232 . t - current time 233 . x - coordinates of the current point 234 . Nc - the number of field components 235 . u - the solution field evaluated at the current point 236 - ctx - a user context, set with `PetscDSSetExactSolution()` or `PetscDSSetExactSolutionTimeDerivative()` 237 238 Level: beginner 239 240 .seealso: `PetscPointFn`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolution()`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolutionTimeDerivative()` 241 S*/ 242 typedef PetscErrorCode PetscPointExactSolutionFn(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx); 243 244 /*S 245 PetscRiemannFn - A prototype of a pointwise function that can be passed to, for example, `PetscDSSetRiemannSolver()` 246 247 Calling Sequence: 248 + dim - the coordinate dimension 249 . Nf - The number of fields 250 . x - The coordinates at a point on the interface 251 . n - The normal vector to the interface 252 . uL - The state vector to the left of the interface 253 . uR - The state vector to the right of the interface 254 . numConstants - number of constant parameters 255 . constants - constant parameters 256 . flux - output array of flux through the interface 257 - ctx - optional user context 258 259 Level: beginner 260 261 .seealso: `PetscPointFn`, `PetscDSSetRiemannSolver()`, `PetscDSGetRiemannSolver()` 262 S*/ 263 typedef void PetscRiemannFn(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx); 264 265 /*S 266 PetscSimplePointFn - A prototype of a simple pointwise function that can be passed to, for example, `DMPlexTransformExtrudeSetNormalFunction()` 267 268 Calling Sequence: 269 + dim - The coordinate dimension of the original mesh (usually a surface) 270 . time - The current time, or 0. 271 . x - The location of the current normal, in the coordinate space of the original mesh 272 . r - The layer number of this point 273 . u - The user provides the computed normal on output 274 - ctx - An optional user context, this context may be obtained by the calling code with `DMGetApplicationContext()` 275 276 Level: beginner 277 278 Developer Note: 279 The handling of `ctx` in the use of such functions may not be ideal since the context is not provided when the function pointer is provided with, for example, `DMSwarmSetCoordinateFunction()` 280 281 .seealso: `PetscPointFn`, `DMPlexTransformExtrudeSetNormalFunction()`, `DMSwarmSetCoordinateFunction()` 282 S*/ 283 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PetscSimplePointFn(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt r, PetscScalar u[], void *ctx); 284 285 PETSC_EXTERN_TYPEDEF typedef PetscSimplePointFn *PetscSimplePointFunc PETSC_DEPRECATED_TYPEDEF(3, 24, 0, "PetscSimplePointFn*", ); 286 PETSC_EXTERN_TYPEDEF typedef PetscPointFn *PetscPointFunc PETSC_DEPRECATED_TYPEDEF(3, 24, 0, "PetscPointFn*", ); 287 PETSC_EXTERN_TYPEDEF typedef PetscPointJacFn *PetscPointJac PETSC_DEPRECATED_TYPEDEF(3, 24, 0, "PetscPointJacFn*", ); 288 PETSC_EXTERN_TYPEDEF typedef PetscBdPointFn *PetscBdPointFunc PETSC_DEPRECATED_TYPEDEF(3, 24, 0, "PetscBdPointFn*", ); 289 PETSC_EXTERN_TYPEDEF typedef PetscBdPointJacFn *PetscBdPointJac PETSC_DEPRECATED_TYPEDEF(3, 24, 0, "PetscBdPointJacFn*", ); 290 PETSC_EXTERN_TYPEDEF typedef PetscRiemannFn *PetscRiemannFunc PETSC_DEPRECATED_TYPEDEF(3, 24, 0, "PetscRiemannFn*", ); 291 292 /*S 293 PetscPointBoundFn - A prototype of a pointwise function that can be passed to, for example, `PetscDSSetLowerBound()` 294 295 Calling Sequence: 296 + dim - the coordinate dimension 297 . t - current time 298 . x - coordinates of the current point 299 . Nc - the number of field components 300 . u - the lower bound evaluated at the current point 301 - ctx - a user context, passed in with, for example, `PetscDSSetLowerBound()` 302 303 Level: beginner 304 305 .seealso: `PetscPointFn`, `PetscDSSetLowerBound()`, `PetscDSSetUpperBound()` 306 S*/ 307 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PetscPointBoundFn(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx); 308