1 #pragma once
2
3 #include <petscds.h>
4 #include <petsc/private/petscimpl.h>
5 #include <petsc/private/hashmap.h>
6
7 PETSC_EXTERN PetscBool PetscDSRegisterAllCalled;
8 PETSC_EXTERN PetscErrorCode PetscDSRegisterAll(void);
9
10 typedef struct _n_DSBoundary *DSBoundary;
11
12 struct _n_DSBoundary {
13 const char *name; /* A unique identifier for the condition */
14 DMLabel label; /* The DMLabel indicating the mesh region over which the condition holds */
15 const char *lname; /* The label name if the label is missing from the DM */
16 PetscInt Nv; /* The number of label values defining the region */
17 PetscInt *values; /* The labels values defining the region */
18 PetscWeakForm wf; /* Holds the pointwise functions defining the form (only for NATURAL conditions) */
19 DMBoundaryConditionType type; /* The type of condition, usually either ESSENTIAL or NATURAL */
20 PetscInt field; /* The field constrained by the condition */
21 PetscInt Nc; /* The number of constrained field components */
22 PetscInt *comps; /* The constrained field components */
23 PetscVoidFn *func; /* Function that provides the boundary values (only for ESSENTIAL conditions) */
24 PetscVoidFn *func_t; /* Function that provides the time derivative of the boundary values (only for ESSENTIAL conditions) */
25 void *ctx; /* The user context for func and func_t */
26 DSBoundary next;
27 };
28
29 typedef struct {
30 PetscCount start; /* Starting entry of the chunk in an array (in bytes) */
31 PetscCount size; /* Current number of entries of the chunk */
32 PetscCount reserved; /* Number of reserved entries in the chunk */
33 } PetscChunk;
34
35 typedef struct {
36 PetscCount size; /* Current number of entries used in array */
37 size_t alloc; /* Number of bytes allocated for array */
38 size_t unitbytes; /* Number of bytes per entry */
39 char *array;
40 } PetscChunkBuffer;
41
42 #define PetscFormKeyHash(key) PetscHashCombine(PetscHashCombine(PetscHashCombine(PetscHashPointer((key).label), PetscHashInt((key).value)), PetscHashInt((key).field)), PetscHashInt((key).part))
43
44 #define PetscFormKeyEqual(k1, k2) (((k1).label == (k2).label) ? ((k1).value == (k2).value) ? ((k1).field == (k2).field) ? ((k1).part == (k2).part) : 0 : 0 : 0)
45
46 static PetscChunk _PetscInvalidChunk = {-1, -1, -1};
47
PETSC_HASH_MAP(HMapForm,PetscFormKey,PetscChunk,PetscFormKeyHash,PetscFormKeyEqual,_PetscInvalidChunk)48 PETSC_HASH_MAP(HMapForm, PetscFormKey, PetscChunk, PetscFormKeyHash, PetscFormKeyEqual, _PetscInvalidChunk)
49
50 /*
51 We sort lexicographically on the structure.
52 Returns
53 -1: left < right
54 0: left = right
55 1: left > right
56 */
57 static inline int Compare_PetscFormKey_Private(const void *left, const void *right, PETSC_UNUSED PetscCtx ctx)
58 {
59 PetscFormKey l = *(const PetscFormKey *)left;
60 PetscFormKey r = *(const PetscFormKey *)right;
61 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)))));
62 }
63
64 typedef struct _PetscWeakFormOps *PetscWeakFormOps;
65 struct _PetscWeakFormOps {
66 PetscErrorCode (*setfromoptions)(PetscWeakForm);
67 PetscErrorCode (*setup)(PetscWeakForm);
68 PetscErrorCode (*view)(PetscWeakForm, PetscViewer);
69 PetscErrorCode (*destroy)(PetscWeakForm);
70 };
71
72 struct _p_PetscWeakForm {
73 PETSCHEADER(struct _PetscWeakFormOps);
74 void *data; /* Implementation object */
75
76 PetscInt Nf; /* The number of fields in the system */
77 PetscChunkBuffer *funcs; /* Buffer holding all function pointers */
78 PetscHMapForm *form; /* Stores function pointers for forms */
79 };
80
81 typedef struct _PetscDSOps *PetscDSOps;
82 struct _PetscDSOps {
83 PetscErrorCode (*setfromoptions)(PetscDS);
84 PetscErrorCode (*setup)(PetscDS);
85 PetscErrorCode (*view)(PetscDS, PetscViewer);
86 PetscErrorCode (*destroy)(PetscDS);
87 };
88
89 struct _p_PetscDS {
90 PETSCHEADER(struct _PetscDSOps);
91 void *data; /* Implementation object */
92 PetscDS *subprobs; /* The subspaces for each dimension */
93 PetscBool setup; /* Flag for setup */
94 PetscInt dimEmbed; /* The real space coordinate dimension */
95 PetscInt Nf; /* The number of solution fields */
96 PetscObject *disc; /* The discretization for each solution field (PetscFE, PetscFV, etc.) */
97 PetscBool *cohesive; /* Flag for cohesive discretization */
98 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 */
99 PetscInt printIntegrate; /* Debugging level for kernels */
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 PetscPointFn **update; // Direct update of field coefficients
110 PetscSimplePointFn **exactSol; // Exact solutions for each field
111 void **exactCtx; // Contexts for the exact solution functions
112 PetscSimplePointFn **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 PetscSimplePointFn **lowerBound; // Lower bounds for each each field
115 void **lowerCtx; // Contexts for the lower bounds functions
116 PetscSimplePointFn **upperBound; // Upper bounds for each each field
117 void **upperCtx; // Contexts for the upper bounds functions
118 PetscInt numConstants; // Number of constants passed to all point functions
119 PetscInt numFuncConstants; // Number of constant passed to an individual point function (like field)
120 PetscScalar *constants; // Array of constants passed to point functions
121 void **ctx; // User contexts for each field
122 /* Computed sizes */
123 PetscInt totDim; /* Total system dimension */
124 PetscInt totComp; /* Total field components */
125 PetscInt *Nc; /* Number of components for each field */
126 PetscInt *Nb; /* Number of basis functions for each field */
127 PetscInt *off; /* Offsets for each field */
128 PetscInt *offDer; /* Derivative offsets for each field */
129 PetscInt *offCohesive[3]; /* Offsets for each field on side s of a cohesive cell */
130 PetscInt *offDerCohesive[3]; /* Derivative offsets for each field on side s of a cohesive cell */
131 PetscTabulation *T; /* Basis function and derivative tabulation for each field */
132 PetscTabulation *Tf; /* Basis function and derivative tabulation for each local face and field */
133 /* Work space */
134 PetscScalar *u; /* Field evaluation */
135 PetscScalar *u_t; /* Field time derivative evaluation */
136 PetscScalar *u_x; /* Field gradient evaluation */
137 PetscScalar *basisReal; /* Workspace for pushforward */
138 PetscScalar *basisDerReal; /* Workspace for derivative pushforward */
139 PetscScalar *testReal; /* Workspace for pushforward */
140 PetscScalar *testDerReal; /* Workspace for derivative pushforward */
141 PetscReal *x; /* Workspace for computing real coordinates */
142 PetscScalar *f0, *f1; /* Point evaluations of weak form residual integrands */
143 PetscScalar *g0, *g1, *g2, *g3; /* Point evaluations of weak form Jacobian integrands */
144 };
145
146 typedef struct {
147 PetscInt dummy; /* */
148 } PetscDS_Basic;
149
150 PETSC_INTERN PetscErrorCode PetscDSGetDiscType_Internal(PetscDS, PetscInt, PetscDiscType *);
151