Actual source code: petscdsimpl.h

  1: #if !defined(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;
 15:   const char *labelname;
 16:   DMBoundaryConditionType type;
 17:   PetscInt    field;
 18:   PetscInt    numcomps;
 19:   PetscInt   *comps;
 20:   void      (*func)(void);
 21:   void      (*func_t)(void);
 22:   PetscInt    numids;
 23:   PetscInt   *ids;
 24:   void       *ctx;
 25:   DSBoundary  next;
 26: };

 28: typedef struct {
 29:   PetscInt start;    /* Starting entry of the chunk in an array (in bytes) */
 30:   PetscInt size;     /* Current number of entries of the chunk */
 31:   PetscInt reserved; /* Number of reserved entries in the chunk */
 32: } PetscChunk;

 34: typedef struct {
 35:   size_t  size;      /* Current number of entries used in array */
 36:   size_t  alloc;     /* Number of bytes allocated for array */
 37:   size_t  unitbytes; /* Number of bytes per entry */
 38:   char   *array;
 39: } PetscChunkBuffer;

 41: #define PetscHashFormKeyHash(key) \
 42:   PetscHashCombine(PetscHashCombine(PetscHashPointer((key).label),PetscHashInt((key).value)),PetscHashInt((key).field))

 44: #define PetscHashFormKeyEqual(k1,k2) \
 45:   (((k1).label == (k2).label) ? \
 46:    ((k1).value == (k2).value) ? \
 47:    ((k1).field == (k2).field) : 0 : 0)

 49: static PetscChunk _PetscInvalidChunk = {-1, -1, -1};

 51: PETSC_HASH_MAP(HMapForm, PetscHashFormKey, PetscChunk, PetscHashFormKeyHash, PetscHashFormKeyEqual, _PetscInvalidChunk)

 53: /*
 54:   We sort lexicographically on the structure.
 55:   Returns
 56:   -1: left < right
 57:    0: left = right
 58:    1: left > right
 59: */
 60: PETSC_STATIC_INLINE int Compare_PetscHashFormKey_Private(const void *left, const void *right, PETSC_UNUSED void *ctx)
 61: {
 62:   PetscHashFormKey l = *(PetscHashFormKey *) left;
 63:   PetscHashFormKey r = *(PetscHashFormKey *) right;
 64:   return (l.label < r.label) ? -1 : ((l.label > r.label) ? 1 :
 65:            ((l.value < r.value) ? -1 : (l.value > r.value) ? 1 :
 66:              ((l.field < r.field) ? -1 : (l.field > r.field))));
 67: }

 69: typedef struct _PetscWeakFormOps *PetscWeakFormOps;
 70: struct _PetscWeakFormOps {
 71:   PetscErrorCode (*setfromoptions)(PetscWeakForm);
 72:   PetscErrorCode (*setup)(PetscWeakForm);
 73:   PetscErrorCode (*view)(PetscWeakForm,PetscViewer);
 74:   PetscErrorCode (*destroy)(PetscWeakForm);
 75: };

 77: struct _p_PetscWeakForm {
 78:   PETSCHEADER(struct _PetscWeakFormOps);
 79:   void *data; /* Implementation object */

 81:   PetscInt          Nf;    /* The number of fields in the system */
 82:   PetscChunkBuffer *funcs; /* Buffer holding all function pointers */
 83:   PetscHMapForm     obj;   /* Scalar integral (like an objective function) */
 84:   PetscHMapForm     f0;    /* Weak form integrands against test function for F */
 85:   PetscHMapForm     f1;    /* Weak form integrands against test function derivative for F */
 86:   PetscHMapForm     g0;    /* Weak form integrands for J = dF/du */
 87:   PetscHMapForm     g1;    /* Weak form integrands for J = dF/du */
 88:   PetscHMapForm     g2;    /* Weak form integrands for J = dF/du */
 89:   PetscHMapForm     g3;    /* Weak form integrands for J = dF/du */
 90:   PetscHMapForm     gp0;   /* Weak form integrands for preconditioner for J */
 91:   PetscHMapForm     gp1;   /* Weak form integrands for preconditioner for J */
 92:   PetscHMapForm     gp2;   /* Weak form integrands for preconditioner for J */
 93:   PetscHMapForm     gp3;   /* Weak form integrands for preconditioner for J */
 94:   PetscHMapForm     gt0;   /* Weak form integrands for dF/du_t */
 95:   PetscHMapForm     gt1;   /* Weak form integrands for dF/du_t */
 96:   PetscHMapForm     gt2;   /* Weak form integrands for dF/du_t */
 97:   PetscHMapForm     gt3;   /* Weak form integrands for dF/du_t */
 98:   PetscHMapForm     bdf0;  /* Weak form boundary integrands F_bd */
 99:   PetscHMapForm     bdf1;  /* Weak form boundary integrands F_bd */
100:   PetscHMapForm     bdg0;  /* Weak form boundary integrands J_bd = dF_bd/du */
101:   PetscHMapForm     bdg1;  /* Weak form boundary integrands J_bd = dF_bd/du */
102:   PetscHMapForm     bdg2;  /* Weak form boundary integrands J_bd = dF_bd/du */
103:   PetscHMapForm     bdg3;  /* Weak form boundary integrands J_bd = dF_bd/du */
104:   PetscHMapForm     bdgp0; /* Weak form integrands for preconditioner for J_bd */
105:   PetscHMapForm     bdgp1; /* Weak form integrands for preconditioner for J_bd */
106:   PetscHMapForm     bdgp2; /* Weak form integrands for preconditioner for J_bd */
107:   PetscHMapForm     bdgp3; /* Weak form integrands for preconditioner for J_bd */
108:   PetscHMapForm     r;     /* Riemann solvers */
109: };

111: typedef struct _PetscDSOps *PetscDSOps;
112: struct _PetscDSOps {
113:   PetscErrorCode (*setfromoptions)(PetscDS);
114:   PetscErrorCode (*setup)(PetscDS);
115:   PetscErrorCode (*view)(PetscDS,PetscViewer);
116:   PetscErrorCode (*destroy)(PetscDS);
117: };

119: struct _p_PetscDS {
120:   PETSCHEADER(struct _PetscDSOps);
121:   void        *data;              /* Implementation object */
122:   PetscDS     *subprobs;          /* The subspaces for each dimension */
123:   PetscBool    setup;             /* Flag for setup */
124:   PetscBool    isHybrid;          /* Flag for hybrid cell (this is crappy, but the only thing I can see to do now) */
125:   PetscInt     dimEmbed;          /* The real space coordinate dimension */
126:   PetscInt     Nf;                /* The number of solution fields */
127:   PetscObject *disc;              /* The discretization for each solution field (PetscFE, PetscFV, etc.) */
128:   /* Equations */
129:   DSBoundary            boundary;      /* Linked list of boundary conditions */
130:   PetscBool             useJacPre;     /* Flag for using the Jacobian preconditioner */
131:   PetscBool            *implicit;      /* Flag for implicit or explicit solve for each field */
132:   PetscInt             *jetDegree;     /* The highest derivative for each field equation, or the k-jet that each discretization needs to tabulate */
133:   PetscWeakForm         wf;            /* The PetscWeakForm holding all pointwise functions */
134:   PetscPointFunc       *update;        /* Direct update of field coefficients */
135:   PetscSimplePointFunc *exactSol;      /* Exact solutions for each field */
136:   void                **exactCtx;      /* Contexts for the exact solution functions */
137:   PetscSimplePointFunc *exactSol_t;    /* Time derivative of the exact solutions for each field */
138:   void                **exactCtx_t;    /* Contexts for the time derivative of the exact solution functions */
139:   PetscInt              numConstants;  /* Number of constants passed to point functions */
140:   PetscScalar          *constants;     /* Array of constants passed to point functions */
141:   void                 **ctx;          /* User contexts for each field */
142:   /* Computed sizes */
143:   PetscInt         totDim;             /* Total system dimension */
144:   PetscInt         totComp;            /* Total field components */
145:   PetscInt        *Nc;                 /* Number of components for each field */
146:   PetscInt        *Nb;                 /* Number of basis functions for each field */
147:   PetscInt        *off;                /* Offsets for each field */
148:   PetscInt        *offDer;             /* Derivative offsets for each field */
149:   PetscTabulation *T;                  /* Basis function and derivative tabulation for each field */
150:   PetscTabulation *Tf;                 /* Basis function and derivative tabulation for each local face and field */
151:   /* Work space */
152:   PetscScalar *u;                      /* Field evaluation */
153:   PetscScalar *u_t;                    /* Field time derivative evaluation */
154:   PetscScalar *u_x;                    /* Field gradient evaluation */
155:   PetscScalar *basisReal;              /* Workspace for pushforward */
156:   PetscScalar *basisDerReal;           /* Workspace for derivative pushforward */
157:   PetscScalar *testReal;               /* Workspace for pushforward */
158:   PetscScalar *testDerReal;            /* Workspace for derivative pushforward */
159:   PetscReal   *x;                      /* Workspace for computing real coordinates */
160:   PetscScalar *f0, *f1;                /* Point evaluations of weak form residual integrands */
161:   PetscScalar *g0, *g1, *g2, *g3;      /* Point evaluations of weak form Jacobian integrands */
162: };

164: typedef struct {
165:   PetscInt dummy; /* */
166: } PetscDS_Basic;

168: PETSC_INTERN PetscErrorCode PetscDSGetDiscType_Internal(PetscDS, PetscInt, PetscDiscType *);

170: #endif