Actual source code: petscdsimpl.h

petsc-main 2021-04-20
Report Typos and Errors
  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;          /* 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 nubmer 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 PetscHashFormKeyHash(key) \
 44:   PetscHashCombine(PetscHashCombine(PetscHashPointer((key).label),PetscHashInt((key).value)),PetscHashInt((key).field))

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

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

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

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

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

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

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

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

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

166: typedef struct {
167:   PetscInt dummy; /* */
168: } PetscDS_Basic;

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

172: #endif