Actual source code: sectionimpl.h
petsc-3.12.0 2019-09-29
1: #if !defined(PETSCSECTIONIMPL_H)
2: #define PETSCSECTIONIMPL_H
4: #include <petscsection.h>
5: #include <petsc/private/petscimpl.h>
7: struct _p_PetscSection {
8: PETSCHEADER(int);
9: PetscInt pStart, pEnd; /* The chart: all points are contained in [pStart, pEnd) */
10: IS perm; /* A permutation of [0, pEnd-pStart) */
11: PetscBool pointMajor; /* True if the offsets are point major, otherwise they are fieldMajor */
12: PetscInt *atlasDof; /* Describes layout of storage, point --> # of values */
13: PetscInt *atlasOff; /* Describes layout of storage, point --> offset into storage */
14: PetscInt maxDof; /* Maximum dof on any point */
15: PetscSection bc; /* Describes constraints, point --> # local dofs which are constrained */
16: PetscInt *bcIndices; /* Local indices for constrained dofs */
17: PetscBool setup;
19: PetscInt numFields; /* The number of fields making up the degrees of freedom */
20: char **fieldNames; /* The field names */
21: PetscInt *numFieldComponents; /* The number of components in each field */
22: PetscSection *field; /* A section describing the layout and constraints for each field */
23: PetscBool useFieldOff; /* Use the field offsets directly for the global section, rather than the point offset */
25: PetscObject clObj; /* Key for the closure (right now we only have one) */
26: PetscSection clSection; /* Section giving the number of points in each closure */
27: IS clPoints; /* Points in each closure */
28: PetscInt clSize; /* The size of a dof closure of a cell, when it is uniform */
29: PetscInt *clPerm; /* A permutation of the cell dof closure, of size clSize */
30: PetscInt *clInvPerm; /* The inverse of clPerm */
31: PetscSectionSym sym; /* Symmetries of the data */
32: };
34: struct _PetscSectionSymOps {
35: PetscErrorCode (*getpoints)(PetscSectionSym,PetscSection,PetscInt,const PetscInt *,const PetscInt **,const PetscScalar **);
36: PetscErrorCode (*destroy)(PetscSectionSym);
37: PetscErrorCode (*view)(PetscSectionSym,PetscViewer);
38: };
40: typedef struct _n_SymWorkLink *SymWorkLink;
42: struct _n_SymWorkLink
43: {
44: SymWorkLink next;
45: const PetscInt **perms;
46: const PetscScalar **rots;
47: PetscInt numPoints;
48: };
50: struct _p_PetscSectionSym {
51: PETSCHEADER(struct _PetscSectionSymOps);
52: void *data;
53: SymWorkLink workin;
54: SymWorkLink workout;
55: };
57: PETSC_EXTERN PetscErrorCode PetscSectionSetClosurePermutation_Internal(PetscSection, PetscObject, PetscInt, PetscCopyMode, PetscInt *);
58: PETSC_EXTERN PetscErrorCode PetscSectionGetClosurePermutation_Internal(PetscSection, PetscObject, PetscInt *, const PetscInt *[]);
59: PETSC_EXTERN PetscErrorCode PetscSectionGetClosureInversePermutation_Internal(PetscSection, PetscObject, PetscInt *, const PetscInt *[]);
60: PETSC_EXTERN PetscErrorCode ISIntersect_Caching_Internal(IS, IS, IS *);
62: #endif