1 /* 2 Objects which encapsulate finite element spaces 3 */ 4 #ifndef PETSCDUALSPACE_H 5 #define PETSCDUALSPACE_H 6 #include <petscdm.h> 7 #include <petscdt.h> 8 #include <petscfetypes.h> 9 #include <petscdstypes.h> 10 #include <petscspace.h> 11 12 /* SUBMANSEC = DUALSPACE */ 13 14 /*S 15 PetscDualSpace - PETSc object that manages the dual space to a linear space, e.g. the space of evaluation functionals at the vertices of a triangle 16 17 Level: beginner 18 19 .seealso: `PetscDualSpaceCreate()`, `PetscSpaceCreate()`, `PetscDualSpaceSetType()`, `PetscDualSpaceType` 20 S*/ 21 typedef struct _p_PetscDualSpace *PetscDualSpace; 22 23 /*MC 24 PetscDualSpaceReferenceCell - The type of reference cell 25 26 Level: beginner 27 28 Note: 29 This is used only for automatic creation of reference cells. A `PetscDualSpace` can accept an arbitrary `DM` for a reference cell. 30 31 .seealso: `PetscSpace` 32 M*/ 33 typedef enum { 34 PETSCDUALSPACE_REFCELL_SIMPLEX, 35 PETSCDUALSPACE_REFCELL_TENSOR 36 } PetscDualSpaceReferenceCell; 37 PETSC_EXTERN const char *const PetscDualSpaceReferenceCells[]; 38 39 /*MC 40 PetscDualSpaceTransformType - The type of function transform 41 42 Level: beginner 43 44 Notes: 45 These transforms, and their inverses, are used to move functions and functionals between the reference element and real space. 46 Suppose that we have a mapping $\phi$ which maps the reference cell to real space, and its Jacobian $J$. If we want to transform function $F$ on the reference element, 47 so that it acts on real space, we use the pushforward transform $\sigma^*$. The pullback $\sigma_*$ is the inverse transform. 48 49 .vb 50 Covariant Piola: $\sigma^*(F) = J^{-T} F \circ \phi^{-1)$ 51 Contravariant Piola: $\sigma^*(F) = 1/|J| J F \circ \phi^{-1)$ 52 .ve 53 54 References: 55 . Rognes, Kirby, and Logg, Efficient Assembly of Hdiv and Hrot Conforming Finite Elements, SISC, 31(6), 4130-4151, arXiv 1205.3085, 2010 56 57 .seealso: `PetscDualSpaceGetDeRahm()` 58 M*/ 59 typedef enum { 60 IDENTITY_TRANSFORM, 61 COVARIANT_PIOLA_TRANSFORM, 62 CONTRAVARIANT_PIOLA_TRANSFORM 63 } PetscDualSpaceTransformType; 64 65 PETSC_EXTERN PetscClassId PETSCDUALSPACE_CLASSID; 66 67 /*J 68 PetscDualSpaceType - String with the name of a PETSc dual space 69 70 Level: beginner 71 72 .seealso: `PetscDualSpaceSetType()`, `PetscDualSpace` 73 J*/ 74 typedef const char *PetscDualSpaceType; 75 #define PETSCDUALSPACELAGRANGE "lagrange" 76 #define PETSCDUALSPACESIMPLE "simple" 77 #define PETSCDUALSPACEREFINED "refined" 78 #define PETSCDUALSPACEBDM "bdm" 79 80 /*MC 81 PETSCDUALSPACEBDM = "bdm" - A `PetscDualSpace` object that encapsulates a dual space for Brezzi-Douglas-Marini elements 82 83 Level: intermediate 84 85 Note: 86 This type is a constructor alias of `PETSCDUALSPACELAGRANGE`. During 87 `PetscDualSpaceSetUp()`, the correct value of `PetscDualSpaceSetFormDegree()` is 88 set for H-div conforming spaces. The type of the dual space is then changed to 89 to `PETSCDUALSPACELAGRANGE`. 90 91 .seealso: `PetscDualSpaceType`, `PetscDualSpaceCreate()`, `PetscDualSpaceSetType()`, `PETSCDUALSPACELAGRANGE`, `PetscDualSpaceSetFormDegree()` 92 M*/ 93 94 PETSC_EXTERN PetscFunctionList PetscDualSpaceList; 95 PETSC_EXTERN PetscErrorCode PetscDualSpaceCreate(MPI_Comm, PetscDualSpace *); 96 PETSC_EXTERN PetscErrorCode PetscDualSpaceDestroy(PetscDualSpace *); 97 PETSC_EXTERN PetscErrorCode PetscDualSpaceDuplicate(PetscDualSpace, PetscDualSpace *); 98 PETSC_EXTERN PetscErrorCode PetscDualSpaceSetType(PetscDualSpace, PetscDualSpaceType); 99 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetType(PetscDualSpace, PetscDualSpaceType *); 100 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetUniform(PetscDualSpace, PetscBool *); 101 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetNumDof(PetscDualSpace, const PetscInt **); 102 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetSection(PetscDualSpace, PetscSection *); 103 PETSC_EXTERN PetscErrorCode PetscDualSpaceSetUp(PetscDualSpace); 104 PETSC_EXTERN PetscErrorCode PetscDualSpaceSetFromOptions(PetscDualSpace); 105 PETSC_EXTERN PetscErrorCode PetscDualSpaceViewFromOptions(PetscDualSpace, PetscObject, const char[]); 106 107 PETSC_EXTERN PetscErrorCode PetscDualSpaceView(PetscDualSpace, PetscViewer); 108 PETSC_EXTERN PetscErrorCode PetscDualSpaceRegister(const char[], PetscErrorCode (*)(PetscDualSpace)); 109 PETSC_EXTERN PetscErrorCode PetscDualSpaceRegisterDestroy(void); 110 111 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetDimension(PetscDualSpace, PetscInt *); 112 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetInteriorDimension(PetscDualSpace, PetscInt *); 113 PETSC_EXTERN PetscErrorCode PetscDualSpaceSetNumComponents(PetscDualSpace, PetscInt); 114 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetNumComponents(PetscDualSpace, PetscInt *); 115 PETSC_EXTERN PetscErrorCode PetscDualSpaceSetOrder(PetscDualSpace, PetscInt); 116 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetOrder(PetscDualSpace, PetscInt *); 117 PETSC_EXTERN PetscErrorCode PetscDualSpaceSetDM(PetscDualSpace, DM); 118 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetDM(PetscDualSpace, DM *); 119 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetFunctional(PetscDualSpace, PetscInt, PetscQuadrature *); 120 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetSymmetries(PetscDualSpace, const PetscInt ****, const PetscScalar ****); 121 122 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetAllData(PetscDualSpace, PetscQuadrature *, Mat *); 123 PETSC_EXTERN PetscErrorCode PetscDualSpaceCreateAllDataDefault(PetscDualSpace, PetscQuadrature *, Mat *); 124 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetInteriorData(PetscDualSpace, PetscQuadrature *, Mat *); 125 PETSC_EXTERN PetscErrorCode PetscDualSpaceCreateInteriorDataDefault(PetscDualSpace, PetscQuadrature *, Mat *); 126 PETSC_EXTERN PetscErrorCode PetscDualSpaceEqual(PetscDualSpace, PetscDualSpace, PetscBool *); 127 128 PETSC_EXTERN PetscErrorCode PetscDualSpaceApplyAll(PetscDualSpace, const PetscScalar *, PetscScalar *); 129 PETSC_EXTERN PetscErrorCode PetscDualSpaceApplyAllDefault(PetscDualSpace, const PetscScalar *, PetscScalar *); 130 PETSC_EXTERN PetscErrorCode PetscDualSpaceApplyInterior(PetscDualSpace, const PetscScalar *, PetscScalar *); 131 PETSC_EXTERN PetscErrorCode PetscDualSpaceApplyInteriorDefault(PetscDualSpace, const PetscScalar *, PetscScalar *); 132 133 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetFormDegree(PetscDualSpace, PetscInt *); 134 PETSC_EXTERN PetscErrorCode PetscDualSpaceSetFormDegree(PetscDualSpace, PetscInt); 135 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetDeRahm(PetscDualSpace, PetscInt *); 136 137 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeGetContinuity(PetscDualSpace, PetscBool *); 138 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeSetContinuity(PetscDualSpace, PetscBool); 139 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeGetTensor(PetscDualSpace, PetscBool *); 140 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeSetTensor(PetscDualSpace, PetscBool); 141 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeGetTrimmed(PetscDualSpace, PetscBool *); 142 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeSetTrimmed(PetscDualSpace, PetscBool); 143 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeGetNodeType(PetscDualSpace, PetscDTNodeType *, PetscBool *, PetscReal *); 144 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeSetNodeType(PetscDualSpace, PetscDTNodeType, PetscBool, PetscReal); 145 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeGetUseMoments(PetscDualSpace, PetscBool *); 146 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeSetUseMoments(PetscDualSpace, PetscBool); 147 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeGetMomentOrder(PetscDualSpace, PetscInt *); 148 PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeSetMomentOrder(PetscDualSpace, PetscInt); 149 150 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetHeightSubspace(PetscDualSpace, PetscInt, PetscDualSpace *); 151 PETSC_EXTERN PetscErrorCode PetscDualSpaceGetPointSubspace(PetscDualSpace, PetscInt, PetscDualSpace *); 152 PETSC_EXTERN PetscErrorCode PetscDualSpaceSimpleSetDimension(PetscDualSpace, PetscInt); 153 PETSC_EXTERN PetscErrorCode PetscDualSpaceSimpleSetFunctional(PetscDualSpace, PetscInt, PetscQuadrature); 154 155 PETSC_EXTERN PetscErrorCode PetscDualSpaceRefinedSetCellSpaces(PetscDualSpace, const PetscDualSpace[]); 156 PETSC_EXTERN PetscErrorCode PetscSpaceCreateSubspace(PetscSpace, PetscDualSpace, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscCopyMode, PetscSpace *); 157 158 #endif 159