1 /* 2 Objects which encapsulate finite element spaces and operations 3 */ 4 #pragma once 5 #include <petscdm.h> 6 #include <petscdt.h> 7 #include <petscfetypes.h> 8 #include <petscdstypes.h> 9 #include <petscspace.h> 10 #include <petscdualspace.h> 11 12 /* SUBMANSEC = FE */ 13 14 /*E 15 PetscFEGeomMode - Describes the type of geometry being encoded. 16 17 Values: 18 + `PETSC_FEGEOM_BASIC` - These are normal dim-cells, with dim == dE, and only bulk data is stored. 19 . `PETSC_FEGEOM_EMBEDDED` - These are dim-cells embedded in a higher dimension, as an embedded manifold, where dim < dE and only bulk data is stored. 20 . `PETSC_FEGEOM_BOUNDARY` - These are dim-cells on the boundary of a dE-mesh, so that dim < dE, and both bulk and s = 1 face data are stored. 21 - `PETSC_FEGEOM_COHESIVE` - These are dim-cells in the interior of a dE-mesh, so that dim < dE, and both bulk and s = 2 face data are stored. 22 23 Level: beginner 24 25 Note: 26 .vb 27 dim - The topological dimension and reference coordinate dimension 28 dE - The real coordinate dimension 29 s - The number of supporting cells for a face 30 .ve 31 32 .seealso: [](ch_dmbase), `PetscFEGeom`, `DM`, `DMPLEX`, `PetscFEGeomCreate()` 33 E*/ 34 typedef enum { 35 PETSC_FEGEOM_BASIC, 36 PETSC_FEGEOM_EMBEDDED, 37 PETSC_FEGEOM_BOUNDARY, 38 PETSC_FEGEOM_COHESIVE 39 } PetscFEGeomMode; 40 41 /*MC 42 PetscFEGeom - Structure for geometric information for `PetscFE` 43 44 Level: intermediate 45 46 Note: 47 This is a struct, not a `PetscObject` 48 49 .seealso: `PetscFE`, `PetscFEGeomCreate()`, `PetscFEGeomDestroy()`, `PetscFEGeomGetChunk()`, `PetscFEGeomRestoreChunk()`, `PetscFEGeomGetPoint()`, `PetscFEGeomGetCellPoint()`, 50 `PetscFEGeomComplete()`, `PetscSpace`, `PetscDualSpace` 51 M*/ 52 typedef struct _n_PetscFEGeom { 53 // We can represent several different types of geometry, which we call modes: 54 // basic: dim == dE, only bulk data 55 // These are normal dim-cells 56 // embedded: dim < dE, only bulk data 57 // These are dim-cells embedded in a higher dimension, as an embedded manifold 58 // boundary: dim < dE, bulk and face data 59 // These are dim-cells on the boundary of a dE-mesh 60 // cohesive: dim < dE, bulk and face data 61 // These are dim-cells in the interior of a dE-mesh 62 // affine: 63 // For all modes, the transforms between real and reference are affine 64 PetscFEGeomMode mode; // The type of geometric data stored 65 PetscBool isAffine; // Flag for affine transforms 66 // Sizes 67 PetscInt dim; // dim: topological dimension and reference coordinate dimension 68 PetscInt dimEmbed; // dE: real coordinate dimension 69 PetscInt numCells; // Nc: Number of mesh points represented in the arrays (points are assumed to be the same DMPolytopeType) 70 PetscInt numPoints; // Np: Number of evaluation points represented in the arrays 71 // Bulk data 72 const PetscReal *xi; // xi[dim] The first point in each cell in reference coordinates 73 PetscReal *v; // v[Nc*Np*dE]: The first point in each cell in real coordinates 74 PetscReal *J; // J[Nc*Np*dE*dE]: The Jacobian of the map from reference to real coordinates (if nonsquare it is completed with orthogonal columns) 75 PetscReal *invJ; // invJ[Nc*Np*dE*dE]: The inverse of the Jacobian of the map from reference to real coordinates (if nonsquare it is completed with orthogonal columns) 76 PetscReal *detJ; // detJ[Nc*Np]: The determinant of J, and if J is non-square it is the volume change 77 // Face data 78 PetscReal *n; // n[Nc*Np*dE]: For faces, the normal to the face in real coordinates, outward for the first supporting cell 79 PetscInt (*face)[4]; // face[Nc][s*2]: For faces, the local face number (cone index) and orientation for this face in each supporting cell 80 PetscReal *suppJ[2]; // sJ[s][Nc*Np*dE*dE]: For faces, the Jacobian for each supporting cell 81 PetscReal *suppInvJ[2]; // sInvJ[s][Nc*Np*dE*dE]: For faces, the inverse Jacobian for each supporting cell 82 PetscReal *suppDetJ[2]; // sdetJ[s][Nc*Np]: For faces, the Jacobian determinant for each supporting cell 83 } PetscFEGeom; 84 85 PETSC_EXTERN PetscErrorCode PetscFEInitializePackage(void); 86 87 PETSC_EXTERN PetscErrorCode PetscFEGeomCreate(PetscQuadrature, PetscInt, PetscInt, PetscFEGeomMode, PetscFEGeom **); 88 PETSC_EXTERN PetscErrorCode PetscFEGeomGetChunk(PetscFEGeom *, PetscInt, PetscInt, PetscFEGeom **); 89 PETSC_EXTERN PetscErrorCode PetscFEGeomRestoreChunk(PetscFEGeom *, PetscInt, PetscInt, PetscFEGeom **); 90 PETSC_EXTERN PetscErrorCode PetscFEGeomGetPoint(PetscFEGeom *, PetscInt, PetscInt, const PetscReal[], PetscFEGeom *); 91 PETSC_EXTERN PetscErrorCode PetscFEGeomGetCellPoint(PetscFEGeom *, PetscInt, PetscInt, PetscFEGeom *); 92 PETSC_EXTERN PetscErrorCode PetscFEGeomComplete(PetscFEGeom *); 93 PETSC_EXTERN PetscErrorCode PetscFEGeomDestroy(PetscFEGeom **); 94 95 PETSC_EXTERN PetscErrorCode PetscDualSpaceApply(PetscDualSpace, PetscInt, PetscReal, PetscFEGeom *, PetscInt, PetscErrorCode (*)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar *, void *), void *, PetscScalar *); 96 PETSC_EXTERN PetscErrorCode PetscDualSpaceApplyDefault(PetscDualSpace, PetscInt, PetscReal, PetscFEGeom *, PetscInt, PetscErrorCode (*)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar *, void *), void *, PetscScalar *); 97 98 PETSC_EXTERN PetscErrorCode PetscDualSpaceTransform(PetscDualSpace, PetscDualSpaceTransformType, PetscBool, PetscFEGeom *, PetscInt, PetscInt, PetscScalar[]); 99 PETSC_EXTERN PetscErrorCode PetscDualSpaceTransformGradient(PetscDualSpace, PetscDualSpaceTransformType, PetscBool, PetscFEGeom *, PetscInt, PetscInt, PetscScalar[]); 100 PETSC_EXTERN PetscErrorCode PetscDualSpaceTransformHessian(PetscDualSpace, PetscDualSpaceTransformType, PetscBool, PetscFEGeom *, PetscInt, PetscInt, PetscScalar[]); 101 PETSC_EXTERN PetscErrorCode PetscDualSpacePullback(PetscDualSpace, PetscFEGeom *, PetscInt, PetscInt, PetscScalar[]); 102 PETSC_EXTERN PetscErrorCode PetscDualSpacePushforward(PetscDualSpace, PetscFEGeom *, PetscInt, PetscInt, PetscScalar[]); 103 PETSC_EXTERN PetscErrorCode PetscDualSpacePushforwardGradient(PetscDualSpace, PetscFEGeom *, PetscInt, PetscInt, PetscScalar[]); 104 PETSC_EXTERN PetscErrorCode PetscDualSpacePushforwardHessian(PetscDualSpace, PetscFEGeom *, PetscInt, PetscInt, PetscScalar[]); 105 106 PETSC_EXTERN PetscClassId PETSCFE_CLASSID; 107 108 /*J 109 PetscFEType - String with the name of a PETSc finite element space 110 111 Level: beginner 112 113 Note: 114 Currently, the classes are concerned with the implementation of element integration 115 116 .seealso: `PetscFESetType()`, `PetscFE` 117 J*/ 118 typedef const char *PetscFEType; 119 #define PETSCFEBASIC "basic" 120 #define PETSCFEOPENCL "opencl" 121 #define PETSCFECOMPOSITE "composite" 122 #define PETSCFEVECTOR "vector" 123 124 PETSC_EXTERN PetscFunctionList PetscFEList; 125 PETSC_EXTERN PetscErrorCode PetscFECreate(MPI_Comm, PetscFE *); 126 PETSC_EXTERN PetscErrorCode PetscFEDestroy(PetscFE *); 127 PETSC_EXTERN PetscErrorCode PetscFESetType(PetscFE, PetscFEType); 128 PETSC_EXTERN PetscErrorCode PetscFEGetType(PetscFE, PetscFEType *); 129 PETSC_EXTERN PetscErrorCode PetscFESetUp(PetscFE); 130 PETSC_EXTERN PetscErrorCode PetscFESetFromOptions(PetscFE); 131 PETSC_EXTERN PetscErrorCode PetscFEViewFromOptions(PetscFE, PetscObject, const char[]); 132 PETSC_EXTERN PetscErrorCode PetscFESetName(PetscFE, const char[]); 133 PETSC_EXTERN PetscErrorCode PetscFECreateVector(PetscFE, PetscInt, PetscBool, PetscBool, PetscFE *); 134 135 PETSC_EXTERN PetscErrorCode PetscFEView(PetscFE, PetscViewer); 136 PETSC_EXTERN PetscErrorCode PetscFERegister(const char[], PetscErrorCode (*)(PetscFE)); 137 PETSC_EXTERN PetscErrorCode PetscFERegisterDestroy(void); 138 PETSC_EXTERN PetscErrorCode PetscFECreateDefault(MPI_Comm, PetscInt, PetscInt, PetscBool, const char[], PetscInt, PetscFE *); 139 PETSC_EXTERN PetscErrorCode PetscFECreateByCell(MPI_Comm, PetscInt, PetscInt, DMPolytopeType, const char[], PetscInt, PetscFE *); 140 PETSC_EXTERN PetscErrorCode PetscFECreateLagrange(MPI_Comm, PetscInt, PetscInt, PetscBool, PetscInt, PetscInt, PetscFE *); 141 PETSC_EXTERN PetscErrorCode PetscFECreateLagrangeByCell(MPI_Comm, PetscInt, PetscInt, DMPolytopeType, PetscInt, PetscInt, PetscFE *); 142 PETSC_EXTERN PetscErrorCode PetscFECreateFromSpaces(PetscSpace, PetscDualSpace, PetscQuadrature, PetscQuadrature, PetscFE *); 143 PETSC_EXTERN PetscErrorCode PetscFELimitDegree(PetscFE, PetscInt, PetscInt, PetscFE *); 144 145 PETSC_EXTERN PetscErrorCode PetscFEGetDimension(PetscFE, PetscInt *); 146 PETSC_EXTERN PetscErrorCode PetscFEGetSpatialDimension(PetscFE, PetscInt *); 147 PETSC_EXTERN PetscErrorCode PetscFESetNumComponents(PetscFE, PetscInt); 148 PETSC_EXTERN PetscErrorCode PetscFEGetNumComponents(PetscFE, PetscInt *); 149 PETSC_EXTERN PetscErrorCode PetscFEGetTileSizes(PetscFE, PetscInt *, PetscInt *, PetscInt *, PetscInt *); 150 PETSC_EXTERN PetscErrorCode PetscFESetTileSizes(PetscFE, PetscInt, PetscInt, PetscInt, PetscInt); 151 PETSC_EXTERN PetscErrorCode PetscFESetBasisSpace(PetscFE, PetscSpace); 152 PETSC_EXTERN PetscErrorCode PetscFEGetBasisSpace(PetscFE, PetscSpace *); 153 PETSC_EXTERN PetscErrorCode PetscFESetDualSpace(PetscFE, PetscDualSpace); 154 PETSC_EXTERN PetscErrorCode PetscFEGetDualSpace(PetscFE, PetscDualSpace *); 155 PETSC_EXTERN PetscErrorCode PetscFESetQuadrature(PetscFE, PetscQuadrature); 156 PETSC_EXTERN PetscErrorCode PetscFEGetQuadrature(PetscFE, PetscQuadrature *); 157 PETSC_EXTERN PetscErrorCode PetscFESetFaceQuadrature(PetscFE, PetscQuadrature); 158 PETSC_EXTERN PetscErrorCode PetscFEGetFaceQuadrature(PetscFE, PetscQuadrature *); 159 PETSC_EXTERN PetscErrorCode PetscFECopyQuadrature(PetscFE, PetscFE); 160 PETSC_EXTERN PetscErrorCode PetscFEGetNumDof(PetscFE, const PetscInt **); 161 162 /* TODO: Need a function to reuse the memory when retabulating the same FE at different points */ 163 PETSC_EXTERN PetscErrorCode PetscFEGetCellTabulation(PetscFE, PetscInt, PetscTabulation *); 164 PETSC_EXTERN PetscErrorCode PetscFEGetFaceTabulation(PetscFE, PetscInt, PetscTabulation *); 165 PETSC_EXTERN PetscErrorCode PetscFEGetFaceCentroidTabulation(PetscFE, PetscTabulation *); 166 PETSC_EXTERN PetscErrorCode PetscFECreateTabulation(PetscFE, PetscInt, PetscInt, const PetscReal[], PetscInt, PetscTabulation *); 167 PETSC_EXTERN PetscErrorCode PetscFEComputeTabulation(PetscFE, PetscInt, const PetscReal[], PetscInt, PetscTabulation); 168 PETSC_EXTERN PetscErrorCode PetscTabulationDestroy(PetscTabulation *); 169 170 PETSC_EXTERN PetscErrorCode PetscFERefine(PetscFE, PetscFE *); 171 PETSC_EXTERN PetscErrorCode PetscFEGetHeightSubspace(PetscFE, PetscInt, PetscFE *); 172 173 PETSC_EXTERN PetscErrorCode PetscFECreateCellGeometry(PetscFE, PetscQuadrature, PetscFEGeom *); 174 PETSC_EXTERN PetscErrorCode PetscFEDestroyCellGeometry(PetscFE, PetscFEGeom *); 175 PETSC_EXTERN PetscErrorCode PetscFEPushforward(PetscFE, PetscFEGeom *, PetscInt, PetscScalar[]); 176 PETSC_EXTERN PetscErrorCode PetscFEPushforwardGradient(PetscFE, PetscFEGeom *, PetscInt, PetscScalar[]); 177 PETSC_EXTERN PetscErrorCode PetscFEPushforwardHessian(PetscFE, PetscFEGeom *, PetscInt, PetscScalar[]); 178 179 PETSC_EXTERN PetscErrorCode PetscFEIntegrate(PetscDS, PetscInt, PetscInt, PetscFEGeom *, const PetscScalar[], PetscDS, const PetscScalar[], PetscScalar[]); 180 PETSC_EXTERN PetscErrorCode PetscFEIntegrateBd(PetscDS, PetscInt, void (*)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt, PetscFEGeom *, const PetscScalar[], PetscDS, const PetscScalar[], PetscScalar[]); 181 PETSC_EXTERN PetscErrorCode PetscFEIntegrateResidual(PetscDS, PetscFormKey, PetscInt, PetscFEGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscScalar[]); 182 PETSC_EXTERN PetscErrorCode PetscFEIntegrateBdResidual(PetscDS, PetscWeakForm, PetscFormKey, PetscInt, PetscFEGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscScalar[]); 183 PETSC_EXTERN PetscErrorCode PetscFEIntegrateHybridResidual(PetscDS, PetscDS, PetscFormKey, PetscInt, PetscInt, PetscFEGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscScalar[]); 184 PETSC_EXTERN PetscErrorCode PetscFEIntegrateJacobian(PetscDS, PetscFEJacobianType, PetscFormKey, PetscInt, PetscFEGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscReal, PetscScalar[]); 185 PETSC_EXTERN PetscErrorCode PetscFEIntegrateBdJacobian(PetscDS, PetscWeakForm, PetscFEJacobianType, PetscFormKey, PetscInt, PetscFEGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscReal, PetscScalar[]); 186 PETSC_EXTERN PetscErrorCode PetscFEIntegrateHybridJacobian(PetscDS, PetscDS, PetscFEJacobianType, PetscFormKey, PetscInt, PetscInt, PetscFEGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscReal, PetscScalar[]); 187 188 PETSC_EXTERN PetscErrorCode PetscFECompositeGetMapping(PetscFE, PetscInt *, const PetscReal *[], const PetscReal *[], const PetscReal *[]); 189 190 PETSC_EXTERN PetscErrorCode PetscFECreateHeightTrace(PetscFE, PetscInt, PetscFE *); 191 PETSC_EXTERN PetscErrorCode PetscFECreatePointTrace(PetscFE, PetscInt, PetscFE *); 192 193 PETSC_EXTERN PetscErrorCode PetscFEOpenCLSetRealType(PetscFE, PetscDataType); 194 PETSC_EXTERN PetscErrorCode PetscFEOpenCLGetRealType(PetscFE, PetscDataType *); 195 196 #ifdef PETSC_HAVE_LIBCEED 197 198 #ifndef PLEXFE_QFUNCTION 199 #define PLEXFE_QFUNCTION(fname, f0_name, f1_name) \ 200 CEED_QFUNCTION(PlexQFunction##fname)(void *ctx, const CeedInt Q, const CeedScalar *const *in, CeedScalar *const *out) \ 201 { \ 202 const CeedScalar *u = in[0], *du = in[1], *qdata = in[2]; \ 203 CeedScalar *v = out[0], *dv = out[1]; \ 204 const PetscInt Nc = 1; \ 205 const PetscInt cdim = 2; \ 206 \ 207 CeedPragmaSIMD for (CeedInt i = 0; i < Q; ++i) \ 208 { \ 209 const PetscInt uOff[2] = {0, Nc}; \ 210 const PetscInt uOff_x[2] = {0, Nc * cdim}; \ 211 const CeedScalar x[2] = {qdata[i + Q * 1], qdata[i + Q * 2]}; \ 212 const CeedScalar invJ[2][2] = { \ 213 {qdata[i + Q * 3], qdata[i + Q * 5]}, \ 214 {qdata[i + Q * 4], qdata[i + Q * 6]} \ 215 }; \ 216 const CeedScalar u_x[2] = {invJ[0][0] * du[i + Q * 0] + invJ[1][0] * du[i + Q * 1], invJ[0][1] * du[i + Q * 0] + invJ[1][1] * du[i + Q * 1]}; \ 217 PetscScalar f0[Nc]; \ 218 PetscScalar f1[Nc * cdim]; \ 219 \ 220 for (PetscInt k = 0; k < Nc; ++k) f0[k] = 0; \ 221 for (PetscInt k = 0; k < Nc * cdim; ++k) f1[k] = 0; \ 222 f0_name(2, 1, 0, uOff, uOff_x, u, NULL, u_x, NULL, NULL, NULL, NULL, NULL, 0.0, x, 0, NULL, f0); \ 223 f1_name(2, 1, 0, uOff, uOff_x, u, NULL, u_x, NULL, NULL, NULL, NULL, NULL, 0.0, x, 0, NULL, f1); \ 224 \ 225 dv[i + Q * 0] = qdata[i + Q * 0] * (invJ[0][0] * f1[0] + invJ[0][1] * f1[1]); \ 226 dv[i + Q * 1] = qdata[i + Q * 0] * (invJ[1][0] * f1[0] + invJ[1][1] * f1[1]); \ 227 v[i] = qdata[i + Q * 0] * f0[0]; \ 228 } \ 229 return CEED_ERROR_SUCCESS; \ 230 } 231 #endif 232 233 #else 234 235 #ifndef PLEXFE_QFUNCTION 236 #define PLEXFE_QFUNCTION(fname, f0_name, f1_name) 237 #endif 238 239 #endif 240