Actual source code: petscdsimpl.h
petsc-3.8.3 2017-12-09
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 PetscErrorCode PetscDSRegisterAll(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: PetscInt numids;
21: PetscInt *ids;
22: void *ctx;
23: DSBoundary next;
24: };
26: typedef struct _PetscDSOps *PetscDSOps;
27: struct _PetscDSOps {
28: PetscErrorCode (*setfromoptions)(PetscDS);
29: PetscErrorCode (*setup)(PetscDS);
30: PetscErrorCode (*view)(PetscDS,PetscViewer);
31: PetscErrorCode (*destroy)(PetscDS);
32: };
34: struct _p_PetscDS {
35: PETSCHEADER(struct _PetscDSOps);
36: void *data; /* Implementation object */
37: PetscBool setup; /* Flag for setup */
38: PetscInt Nf; /* The number of solution fields */
39: PetscBool *implicit; /* Flag for implicit or explicit solve */
40: PetscBool *adjacency; /* Flag for variable influence */
41: PetscObject *disc; /* The discretization for each solution field (PetscFE, PetscFV, etc.) */
42: PetscPointFunc *obj; /* Scalar integral (like an objective function) */
43: PetscPointFunc *f; /* Weak form integrands for F, f_0, f_1 */
44: PetscPointJac *g; /* Weak form integrands for J = dF/du, g_0, g_1, g_2, g_3 */
45: PetscPointJac *gp; /* Weak form integrands for preconditioner for J, g_0, g_1, g_2, g_3 */
46: PetscPointJac *gt; /* Weak form integrands for dF/du_t, g_0, g_1, g_2, g_3 */
47: PetscBdPointFunc *fBd; /* Weak form boundary integrands F_bd, f_0, f_1 */
48: PetscBdPointJac *gBd; /* Weak form boundary integrands J_bd = dF_bd/du, g_0, g_1, g_2, g_3 */
49: PetscRiemannFunc *r; /* Riemann solvers */
50: PetscPointFunc *update; /* Direct update of field coefficients */
51: PetscSimplePointFunc *exactSol; /* Exact solutions for each field */
52: PetscInt numConstants; /* Number of constants passed to point functions */
53: PetscScalar *constants; /* Array of constants passed to point functions */
54: void **ctx; /* User contexts for each field */
55: PetscInt dim; /* The spatial dimension */
56: /* Computed sizes */
57: PetscInt totDim; /* Total system dimension */
58: PetscInt totComp; /* Total field components */
59: PetscInt *Nc; /* Number of components for each field */
60: PetscInt *Nb; /* Number of basis functions for each field */
61: PetscInt *off; /* Offsets for each field */
62: PetscInt *offDer; /* Derivative offsets for each field */
63: PetscReal **basis; /* Default basis tabulation for each field */
64: PetscReal **basisDer; /* Default basis derivative tabulation for each field */
65: PetscReal **basisFace; /* Basis tabulation for each local face and field */
66: PetscReal **basisDerFace; /* Basis derivative tabulation for each local face and field */
67: /* Work space */
68: PetscScalar *u; /* Field evaluation */
69: PetscScalar *u_t; /* Field time derivative evaluation */
70: PetscScalar *u_x; /* Field gradient evaluation */
71: PetscScalar *refSpaceDer; /* Workspace for computing derivative in the reference coordinates */
72: PetscReal *x; /* Workspace for computing real coordinates */
73: PetscScalar *f0, *f1; /* Point evaluations of weak form residual integrands */
74: PetscScalar *g0, *g1, *g2, *g3; /* Point evaluations of weak form Jacobian integrands */
75: DSBoundary boundary; /* Linked list of boundary conditions */
76: };
78: typedef struct {
79: PetscInt dummy; /* */
80: } PetscDS_Basic;
82: #endif