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