1 #pragma once
2
3 #include <petscsection.h>
4 #include <petsc/private/petscimpl.h>
5 #include <petsc/private/hashmap.h>
6
7 typedef struct PetscSectionClosurePermKey {
8 PetscInt depth, size;
9 } PetscSectionClosurePermKey;
10
11 typedef struct {
12 PetscInt *perm, *invPerm;
13 } PetscSectionClosurePermVal;
14
PetscSectionClosurePermHash(PetscSectionClosurePermKey k)15 static inline PetscHash_t PetscSectionClosurePermHash(PetscSectionClosurePermKey k)
16 {
17 return PetscHashCombine(PetscHashInt(k.depth), PetscHashInt(k.size));
18 }
19
PetscSectionClosurePermEqual(PetscSectionClosurePermKey k1,PetscSectionClosurePermKey k2)20 static inline int PetscSectionClosurePermEqual(PetscSectionClosurePermKey k1, PetscSectionClosurePermKey k2)
21 {
22 return k1.depth == k2.depth && k1.size == k2.size;
23 }
24
25 static PetscSectionClosurePermVal PetscSectionClosurePermVal_Empty = {PETSC_NULLPTR, PETSC_NULLPTR};
26 PETSC_HASH_MAP(ClPerm, PetscSectionClosurePermKey, PetscSectionClosurePermVal, PetscSectionClosurePermHash, PetscSectionClosurePermEqual, PetscSectionClosurePermVal_Empty)
27
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), so perm[i] is the ith permuted point */
32 PetscBT blockStarts; /* If present, bit is high for each point that starts a block */
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;
41
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 */
48
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 };
55
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 };
63
64 typedef struct _n_SymWorkLink *SymWorkLink;
65
66 struct _n_SymWorkLink {
67 SymWorkLink next;
68 const PetscInt **perms;
69 const PetscScalar **rots;
70 PetscInt numPoints;
71 };
72
73 struct _p_PetscSectionSym {
74 PETSCHEADER(struct _PetscSectionSymOps);
75 void *data;
76 SymWorkLink workin;
77 SymWorkLink workout;
78 };
79
80 PETSC_SINGLE_LIBRARY_INTERN PetscErrorCode PetscSectionCopy_Internal(PetscSection, PetscSection, PetscBT);
81 PETSC_SINGLE_LIBRARY_INTERN PetscErrorCode PetscSectionSetClosurePermutation_Internal(PetscSection, PetscObject, PetscInt, PetscInt, PetscCopyMode, PetscInt *);
82 PETSC_SINGLE_LIBRARY_INTERN PetscErrorCode PetscSectionGetClosureInversePermutation_Internal(PetscSection, PetscObject, PetscInt, PetscInt, const PetscInt *[]);
83 PETSC_SINGLE_LIBRARY_INTERN PetscErrorCode ISIntersect_Caching_Internal(IS, IS, IS *);
84 #if defined(PETSC_HAVE_HDF5)
85 PETSC_INTERN PetscErrorCode PetscSectionView_HDF5_Internal(PetscSection, PetscViewer);
86 PETSC_INTERN PetscErrorCode PetscSectionLoad_HDF5_Internal(PetscSection, PetscViewer);
87 #endif
88 PETSC_INTERN PetscErrorCode PetscSectionArrayView_ASCII_Internal(PetscSection, void *, PetscDataType, PetscViewer);
89
PetscSectionCheckConstraints_Private(PetscSection s)90 static inline PetscErrorCode PetscSectionCheckConstraints_Private(PetscSection s)
91 {
92 PetscFunctionBegin;
93 if (!s->bc) {
94 PetscCall(PetscSectionCreate(PETSC_COMM_SELF, &s->bc));
95 PetscCall(PetscSectionSetChart(s->bc, s->pStart, s->pEnd));
96 }
97 PetscFunctionReturn(PETSC_SUCCESS);
98 }
99
100 // Call this if you directly modify atlasDof so that maxDof gets recalculated on next PetscSectionGetMaxDof()
PetscSectionInvalidateMaxDof_Internal(PetscSection s)101 static inline PetscErrorCode PetscSectionInvalidateMaxDof_Internal(PetscSection s)
102 {
103 s->maxDof = PETSC_INT_MIN;
104 return PETSC_SUCCESS;
105 }
106
107 #if defined(PETSC_CLANG_STATIC_ANALYZER)
108 void PetscSectionCheckValidField(PetscInt, PetscInt);
109 void PetscSectionCheckValidFieldComponent(PetscInt, PetscInt);
110 #else
111 #define PetscSectionCheckValid_(description, item, nitems) \
112 do { \
113 PetscCheck(((item) >= 0) && ((item) < (nitems)), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Invalid " description " %" PetscInt_FMT "; not in [0, %" PetscInt_FMT ")", (item), (nitems)); \
114 } while (0)
115
116 #define PetscSectionCheckValidFieldComponent(comp, nfieldcomp) PetscSectionCheckValid_("section field component", comp, nfieldcomp)
117
118 #define PetscSectionCheckValidField(field, nfields) PetscSectionCheckValid_("field number", field, nfields)
119 #endif /* PETSC_CLANG_STATIC_ANALYZER */
120