Actual source code: sectionimpl.h
1: #if !defined(PETSCSECTIONIMPL_H)
2: #define PETSCSECTIONIMPL_H
4: #include <petscsection.h>
5: #include <petsc/private/petscimpl.h>
6: #include <petsc/private/hashmap.h>
8: typedef struct PetscSectionClosurePermKey {
9: PetscInt depth, size;
10: } PetscSectionClosurePermKey;
12: typedef struct {
13: PetscInt *perm, *invPerm;
14: } PetscSectionClosurePermVal;
16: static inline PetscHash_t PetscSectionClosurePermHash(PetscSectionClosurePermKey k)
17: {
18: return PetscHashCombine(PetscHashInt(k.depth), PetscHashInt(k.size));
19: }
21: static inline int PetscSectionClosurePermEqual(PetscSectionClosurePermKey k1, PetscSectionClosurePermKey k2)
22: {
23: return k1.depth == k2.depth && k1.size == k2.size;
24: }
26: static PetscSectionClosurePermVal PetscSectionClosurePermVal_Empty = {NULL, NULL};
27: PETSC_HASH_MAP(ClPerm, PetscSectionClosurePermKey, PetscSectionClosurePermVal, PetscSectionClosurePermHash, PetscSectionClosurePermEqual, PetscSectionClosurePermVal_Empty)
29: struct _p_PetscSection {
30: PETSCHEADER(int);
31: PetscInt pStart, pEnd; /* The chart: all points are contained in [pStart, pEnd) */
32: IS perm; /* A permutation of [0, pEnd-pStart) */
33: PetscBool pointMajor; /* True if the offsets are point major, otherwise they are fieldMajor */
34: PetscBool includesConstraints; /* True if constrained dofs are included when computing offsets */
35: PetscInt *atlasDof; /* Describes layout of storage, point --> # of values */
36: PetscInt *atlasOff; /* Describes layout of storage, point --> offset into storage */
37: PetscInt maxDof; /* Maximum dof on any point */
38: PetscSection bc; /* Describes constraints, point --> # local dofs which are constrained */
39: PetscInt *bcIndices; /* Local indices for constrained dofs */
40: PetscBool setup;
42: PetscInt numFields; /* The number of fields making up the degrees of freedom */
43: char **fieldNames; /* The field names */
44: PetscInt *numFieldComponents; /* The number of components in each field */
45: PetscSection *field; /* A section describing the layout and constraints for each field */
46: PetscBool useFieldOff; /* Use the field offsets directly for the global section, rather than the point offset */
47: char ***compNames; /* The component names */
49: PetscObject clObj; /* Key for the closure (right now we only have one) */
50: PetscClPerm clHash; /* Hash of (depth, size) to perm and invPerm */
51: PetscSection clSection; /* Section giving the number of points in each closure */
52: IS clPoints; /* Points in each closure */
53: PetscSectionSym sym; /* Symmetries of the data */
54: };
56: struct _PetscSectionSymOps {
57: PetscErrorCode (*getpoints)(PetscSectionSym,PetscSection,PetscInt,const PetscInt *,const PetscInt **,const PetscScalar **);
58: PetscErrorCode (*distribute)(PetscSectionSym,PetscSF,PetscSectionSym*);
59: PetscErrorCode (*copy)(PetscSectionSym,PetscSectionSym);
60: PetscErrorCode (*destroy)(PetscSectionSym);
61: PetscErrorCode (*view)(PetscSectionSym,PetscViewer);
62: };
64: typedef struct _n_SymWorkLink *SymWorkLink;
66: struct _n_SymWorkLink
67: {
68: SymWorkLink next;
69: const PetscInt **perms;
70: const PetscScalar **rots;
71: PetscInt numPoints;
72: };
74: struct _p_PetscSectionSym {
75: PETSCHEADER(struct _PetscSectionSymOps);
76: void *data;
77: SymWorkLink workin;
78: SymWorkLink workout;
79: };
81: PETSC_EXTERN PetscErrorCode PetscSectionSetClosurePermutation_Internal(PetscSection, PetscObject, PetscInt, PetscInt, PetscCopyMode, PetscInt *);
82: PETSC_EXTERN PetscErrorCode PetscSectionGetClosurePermutation_Internal(PetscSection, PetscObject, PetscInt, PetscInt, const PetscInt *[]);
83: PETSC_EXTERN PetscErrorCode PetscSectionGetClosureInversePermutation_Internal(PetscSection, PetscObject, PetscInt, PetscInt, const PetscInt *[]);
84: PETSC_EXTERN PetscErrorCode ISIntersect_Caching_Internal(IS, IS, IS *);
85: #if defined(PETSC_HAVE_HDF5)
86: PETSC_INTERN PetscErrorCode PetscSectionView_HDF5_Internal(PetscSection, PetscViewer);
87: PETSC_INTERN PetscErrorCode PetscSectionLoad_HDF5_Internal(PetscSection, PetscViewer);
88: #endif
90: static inline PetscErrorCode PetscSectionCheckConstraints_Private(PetscSection s)
91: {
92: if (!s->bc) {
93: PetscSectionCreate(PETSC_COMM_SELF, &s->bc);
94: PetscSectionSetChart(s->bc, s->pStart, s->pEnd);
95: }
96: return 0;
97: }
99: /* Call this if you directly modify atlasDof so that maxDof gets recalculated on next PetscSectionGetMaxDof() */
100: static inline PetscErrorCode PetscSectionInvalidateMaxDof_Internal(PetscSection s)
101: {
102: s->maxDof = PETSC_MIN_INT;
103: return 0;
104: }
106: #if defined(PETSC_CLANG_STATIC_ANALYZER)
107: # define PetscSectionCheckValidField(x,y)
108: # define PetscSectionCheckValidFieldComponent(x,y)
109: #else
110: # define PetscSectionCheckValid_(description,item,nitems) do { \
112: } while (0)
114: # define PetscSectionCheckValidFieldComponent(comp,nfieldcomp) PetscSectionCheckValid_("section field component",comp,nfieldcomp)
116: # define PetscSectionCheckValidField(field,nfields) PetscSectionCheckValid_("field number",field,nfields)
118: #endif /* PETSC_CLANG_STATIC_ANALYZER */
120: #endif /* PETSCSECTIONIMPL_H */