1: #if !defined(PETSCDSIMPL_H)
2: #define PETSCDSIMPL_H 4: #include <petscds.h> 5: #include <petsc/private/petscimpl.h> 7: PETSC_EXTERN PetscBool PetscDSRegisterAllCalled;
8: PETSC_EXTERN PetscErrorCodePetscDSRegisterAll(void);
10: typedef struct _n_DSBoundary *DSBoundary;
12: struct _n_DSBoundary {
13: const char *name;
14: const char *labelname;
15: DMBoundaryConditionType type;
16: PetscInt field;
17: PetscInt numcomps;
18: PetscInt *comps;
19: void (*func)(void);
20: void (*func_t)(void);
21: PetscInt numids;
22: PetscInt *ids;
23: void *ctx;
24: DSBoundary next;
25: };
27: typedef struct _PetscDSOps *PetscDSOps;
28: struct _PetscDSOps {
29: PetscErrorCode (*setfromoptions)(PetscDS);
30: PetscErrorCode (*setup)(PetscDS);
31: PetscErrorCode (*view)(PetscDS,PetscViewer);
32: PetscErrorCode (*destroy)(PetscDS);
33: };
35: struct _p_PetscDS {
36: PETSCHEADER(struct _PetscDSOps);
37: void *data; /* Implementation object */
38: PetscDS *subprobs; /* The subspaces for each dimension */
39: PetscBool setup; /* Flag for setup */
40: PetscBool isHybrid; /* Flag for hybrid cell (this is crappy, but the only thing I can see to do now) */
41: PetscInt dimEmbed; /* The real space coordinate dimension */
42: PetscInt Nf; /* The number of solution fields */
43: PetscObject *disc; /* The discretization for each solution field (PetscFE, PetscFV, etc.) */
44: /* Equations */
45: DSBoundary boundary; /* Linked list of boundary conditions */
46: PetscBool useJacPre; /* Flag for using the Jacobian preconditioner */
47: PetscBool *implicit; /* Flag for implicit or explicit solve for each field */
48: PetscInt *jetDegree; /* The highest derivative for each field equation, or the k-jet that each discretization needs to tabulate */
49: PetscPointFunc *obj; /* Scalar integral (like an objective function) */
50: PetscPointFunc *f; /* Weak form integrands for F, f_0, f_1 */
51: PetscPointJac *g; /* Weak form integrands for J = dF/du, g_0, g_1, g_2, g_3 */
52: PetscPointJac *gp; /* Weak form integrands for preconditioner for J, g_0, g_1, g_2, g_3 */
53: PetscPointJac *gt; /* Weak form integrands for dF/du_t, g_0, g_1, g_2, g_3 */
54: PetscBdPointFunc *fBd; /* Weak form boundary integrands F_bd, f_0, f_1 */
55: PetscBdPointJac *gBd; /* Weak form boundary integrands J_bd = dF_bd/du, g_0, g_1, g_2, g_3 */
56: PetscBdPointJac *gpBd; /* Weak form integrands for preconditioner for J_bd, g_0, g_1, g_2, g_3 */
57: PetscRiemannFunc *r; /* Riemann solvers */
58: PetscPointFunc *update; /* Direct update of field coefficients */
59: PetscSimplePointFunc *exactSol; /* Exact solutions for each field */
60: void **exactCtx; /* Contexts for the exact solution functions */
61: PetscSimplePointFunc *exactSol_t; /* Time derivative of the exact solutions for each field */
62: void **exactCtx_t; /* Contexts for the time derivative of the exact solution functions */
63: PetscInt numConstants; /* Number of constants passed to point functions */
64: PetscScalar *constants; /* Array of constants passed to point functions */
65: void **ctx; /* User contexts for each field */
66: /* Computed sizes */
67: PetscInt totDim; /* Total system dimension */
68: PetscInt totComp; /* Total field components */
69: PetscInt *Nc; /* Number of components for each field */
70: PetscInt *Nb; /* Number of basis functions for each field */
71: PetscInt *off; /* Offsets for each field */
72: PetscInt *offDer; /* Derivative offsets for each field */
73: PetscTabulation *T; /* Basis function and derivative tabulation for each field */
74: PetscTabulation *Tf; /* Basis function and derivative tabulation for each local face and field */
75: /* Work space */
76: PetscScalar *u; /* Field evaluation */
77: PetscScalar *u_t; /* Field time derivative evaluation */
78: PetscScalar *u_x; /* Field gradient evaluation */
79: PetscScalar *basisReal; /* Workspace for pushforward */
80: PetscScalar *basisDerReal; /* Workspace for derivative pushforward */
81: PetscScalar *testReal; /* Workspace for pushforward */
82: PetscScalar *testDerReal; /* Workspace for derivative pushforward */
83: PetscReal *x; /* Workspace for computing real coordinates */
84: PetscScalar *f0, *f1; /* Point evaluations of weak form residual integrands */
85: PetscScalar *g0, *g1, *g2, *g3; /* Point evaluations of weak form Jacobian integrands */
86: };
88: typedef struct {
89: PetscInt dummy; /* */
90: } PetscDS_Basic;
92: PETSC_INTERN PetscErrorCode PetscDSGetDiscType_Internal(PetscDS, PetscInt, PetscDiscType *);
94: #endif