Actual source code: petscdsimpl.h
1: #ifndef PETSCDSIMPL_H
2: #define PETSCDSIMPL_H
4: #include <petscds.h>
5: #include <petsc/private/petscimpl.h>
6: #include <petsc/private/hashmap.h>
8: PETSC_EXTERN PetscBool PetscDSRegisterAllCalled;
9: PETSC_EXTERN PetscErrorCode PetscDSRegisterAll(void);
11: typedef struct _n_DSBoundary *DSBoundary;
13: struct _n_DSBoundary {
14: const char *name; /* A unique identifier for the condition */
15: DMLabel label; /* The DMLabel indicating the mesh region over which the condition holds */
16: const char *lname; /* The label name if the label is missing from the DM */
17: PetscInt Nv; /* The number of label values defining the region */
18: PetscInt *values; /* The labels values defining the region */
19: PetscWeakForm wf; /* Holds the pointwise functions defining the form (only for NATURAL conditions) */
20: DMBoundaryConditionType type; /* The type of condition, usually either ESSENTIAL or NATURAL */
21: PetscInt field; /* The field constrained by the condition */
22: PetscInt Nc; /* The number of constrained field components */
23: PetscInt *comps; /* The constrained field components */
24: void (*func)(void); /* Function that provides the boundary values (only for ESSENTIAL conditions) */
25: void (*func_t)(void); /* Function that provides the time derivative of the boundary values (only for ESSENTIAL conditions) */
26: void *ctx; /* The user context for func and func_t */
27: DSBoundary next;
28: };
30: typedef struct {
31: PetscInt start; /* Starting entry of the chunk in an array (in bytes) */
32: PetscInt size; /* Current number of entries of the chunk */
33: PetscInt reserved; /* Number of reserved entries in the chunk */
34: } PetscChunk;
36: typedef struct {
37: size_t size; /* Current number of entries used in array */
38: size_t alloc; /* Number of bytes allocated for array */
39: size_t unitbytes; /* Number of bytes per entry */
40: char *array;
41: } PetscChunkBuffer;
43: #define PetscFormKeyHash(key) PetscHashCombine(PetscHashCombine(PetscHashCombine(PetscHashPointer((key).label), PetscHashInt((key).value)), PetscHashInt((key).field)), PetscHashInt((key).part))
45: #define PetscFormKeyEqual(k1, k2) (((k1).label == (k2).label) ? ((k1).value == (k2).value) ? ((k1).field == (k2).field) ? ((k1).part == (k2).part) : 0 : 0 : 0)
47: static PetscChunk _PetscInvalidChunk = {-1, -1, -1};
49: PETSC_HASH_MAP(HMapForm, PetscFormKey, PetscChunk, PetscFormKeyHash, PetscFormKeyEqual, _PetscInvalidChunk)
51: /*
52: We sort lexicographically on the structure.
53: Returns
54: -1: left < right
55: 0: left = right
56: 1: left > right
57: */
58: static inline int Compare_PetscFormKey_Private(const void *left, const void *right, PETSC_UNUSED void *ctx)
59: {
60: PetscFormKey l = *(const PetscFormKey *)left;
61: PetscFormKey r = *(const PetscFormKey *)right;
62: return (l.label < r.label) ? -1 : ((l.label > r.label) ? 1 : ((l.value < r.value) ? -1 : (l.value > r.value) ? 1 : ((l.field < r.field) ? -1 : (l.field > r.field) ? 1 : ((l.part < r.part) ? -1 : (l.part > r.part)))));
63: }
65: typedef struct _PetscWeakFormOps *PetscWeakFormOps;
66: struct _PetscWeakFormOps {
67: PetscErrorCode (*setfromoptions)(PetscWeakForm);
68: PetscErrorCode (*setup)(PetscWeakForm);
69: PetscErrorCode (*view)(PetscWeakForm, PetscViewer);
70: PetscErrorCode (*destroy)(PetscWeakForm);
71: };
73: struct _p_PetscWeakForm {
74: PETSCHEADER(struct _PetscWeakFormOps);
75: void *data; /* Implementation object */
77: PetscInt Nf; /* The number of fields in the system */
78: PetscChunkBuffer *funcs; /* Buffer holding all function pointers */
79: PetscHMapForm *form; /* Stores function pointers for forms */
80: };
82: typedef struct _PetscDSOps *PetscDSOps;
83: struct _PetscDSOps {
84: PetscErrorCode (*setfromoptions)(PetscDS);
85: PetscErrorCode (*setup)(PetscDS);
86: PetscErrorCode (*view)(PetscDS, PetscViewer);
87: PetscErrorCode (*destroy)(PetscDS);
88: };
90: struct _p_PetscDS {
91: PETSCHEADER(struct _PetscDSOps);
92: void *data; /* Implementation object */
93: PetscDS *subprobs; /* The subspaces for each dimension */
94: PetscBool setup; /* Flag for setup */
95: PetscInt dimEmbed; /* The real space coordinate dimension */
96: PetscInt Nf; /* The number of solution fields */
97: PetscObject *disc; /* The discretization for each solution field (PetscFE, PetscFV, etc.) */
98: PetscBool *cohesive; /* Flag for cohesive discretization */
99: PetscBool isCohesive; /* We are on a cohesive cell, meaning lower dimensional FE used on a 0-volume cell. Normal fields appear on both endcaps, whereas cohesive field only appear once in the middle */
100: /* Quadrature */
101: PetscBool forceQuad; /* Flag to force matching quadratures in discretizations */
102: IS *quadPerm[DM_NUM_POLYTOPES]; /* qP[ct][o]: q point permutation for orientation o of integ domain */
103: /* Equations */
104: DSBoundary boundary; /* Linked list of boundary conditions */
105: PetscBool useJacPre; /* Flag for using the Jacobian preconditioner */
106: PetscBool *implicit; /* Flag for implicit or explicit solve for each field */
107: PetscInt *jetDegree; /* The highest derivative for each field equation, or the k-jet that each discretization needs to tabulate */
108: PetscWeakForm wf; /* The PetscWeakForm holding all pointwise functions */
109: PetscPointFunc *update; /* Direct update of field coefficients */
110: PetscSimplePointFunc *exactSol; /* Exact solutions for each field */
111: void **exactCtx; /* Contexts for the exact solution functions */
112: PetscSimplePointFunc *exactSol_t; /* Time derivative of the exact solutions for each field */
113: void **exactCtx_t; /* Contexts for the time derivative of the exact solution functions */
114: PetscInt numConstants; /* Number of constants passed to point functions */
115: PetscScalar *constants; /* Array of constants passed to point functions */
116: void **ctx; /* User contexts for each field */
117: /* Computed sizes */
118: PetscInt totDim; /* Total system dimension */
119: PetscInt totComp; /* Total field components */
120: PetscInt *Nc; /* Number of components for each field */
121: PetscInt *Nb; /* Number of basis functions for each field */
122: PetscInt *off; /* Offsets for each field */
123: PetscInt *offDer; /* Derivative offsets for each field */
124: PetscInt *offCohesive[3]; /* Offsets for each field on side s of a cohesive cell */
125: PetscInt *offDerCohesive[3]; /* Derivative offsets for each field on side s of a cohesive cell */
126: PetscTabulation *T; /* Basis function and derivative tabulation for each field */
127: PetscTabulation *Tf; /* Basis function and derivative tabulation for each local face and field */
128: /* Work space */
129: PetscScalar *u; /* Field evaluation */
130: PetscScalar *u_t; /* Field time derivative evaluation */
131: PetscScalar *u_x; /* Field gradient evaluation */
132: PetscScalar *basisReal; /* Workspace for pushforward */
133: PetscScalar *basisDerReal; /* Workspace for derivative pushforward */
134: PetscScalar *testReal; /* Workspace for pushforward */
135: PetscScalar *testDerReal; /* Workspace for derivative pushforward */
136: PetscReal *x; /* Workspace for computing real coordinates */
137: PetscScalar *f0, *f1; /* Point evaluations of weak form residual integrands */
138: PetscScalar *g0, *g1, *g2, *g3; /* Point evaluations of weak form Jacobian integrands */
139: };
141: typedef struct {
142: PetscInt dummy; /* */
143: } PetscDS_Basic;
145: PETSC_INTERN PetscErrorCode PetscDSGetDiscType_Internal(PetscDS, PetscInt, PetscDiscType *);
147: #endif