Actual source code: petscfe.h
petsc-3.9.4 2018-09-11
1: /*
2: Objects which encapsulate finite element spaces and operations
3: */
6: #include <petscdm.h>
7: #include <petscdt.h>
8: #include <petscfetypes.h>
9: #include <petscdstypes.h>
11: /* Assuming dim <= 3 */
12: typedef struct {
13: PetscReal v0[3];
14: PetscReal J[9];
15: PetscReal invJ[9];
16: PetscReal detJ;
17: PetscInt dim;
18: PetscInt dimEmbed;
19: } PetscFECellGeom;
21: typedef struct {
22: PetscReal v0[3];
23: PetscReal J[9]; /* Size could be reduced to 6 */
24: PetscReal invJ[2][9]; /* invJ in the adjacent cells */
25: PetscReal detJ; /* Not really the determinant of J, but \sqrt{\det{J^T J}} */
26: PetscReal n[3];
27: PetscInt dim;
28: PetscInt dimEmbed;
29: PetscInt face[2]; /* The local face numbers in the adjacent cells */
30: } PetscFEFaceGeom;
32: PETSC_EXTERN PetscErrorCode PetscFEInitializePackage(void);
34: PETSC_EXTERN PetscClassId PETSCSPACE_CLASSID;
36: /*J
37: PetscSpaceType - String with the name of a PETSc linear space
39: Level: beginner
41: .seealso: PetscSpaceSetType(), PetscSpace
42: J*/
43: typedef const char *PetscSpaceType;
44: #define PETSCSPACEPOLYNOMIAL "poly"
45: #define PETSCSPACEPOINT "point"
47: PETSC_EXTERN PetscFunctionList PetscSpaceList;
48: PETSC_EXTERN PetscErrorCode PetscSpaceCreate(MPI_Comm, PetscSpace *);
49: PETSC_EXTERN PetscErrorCode PetscSpaceDestroy(PetscSpace *);
50: PETSC_EXTERN PetscErrorCode PetscSpaceSetType(PetscSpace, PetscSpaceType);
51: PETSC_EXTERN PetscErrorCode PetscSpaceGetType(PetscSpace, PetscSpaceType *);
52: PETSC_EXTERN PetscErrorCode PetscSpaceSetUp(PetscSpace);
53: PETSC_EXTERN PetscErrorCode PetscSpaceSetFromOptions(PetscSpace);
54: PETSC_STATIC_INLINE PetscErrorCode PetscSpaceViewFromOptions(PetscSpace A,PetscObject B,const char name[]) {return PetscObjectViewFromOptions((PetscObject)A,B,name);}
56: PETSC_EXTERN PetscErrorCode PetscSpaceView(PetscSpace,PetscViewer);
57: PETSC_EXTERN PetscErrorCode PetscSpaceRegister(const char [], PetscErrorCode (*)(PetscSpace));
58: PETSC_EXTERN PetscErrorCode PetscSpaceRegisterDestroy(void);
60: PETSC_EXTERN PetscErrorCode PetscSpaceGetDimension(PetscSpace, PetscInt *);
61: PETSC_EXTERN PetscErrorCode PetscSpaceSetNumComponents(PetscSpace, PetscInt);
62: PETSC_EXTERN PetscErrorCode PetscSpaceGetNumComponents(PetscSpace, PetscInt *);
63: PETSC_EXTERN PetscErrorCode PetscSpaceSetOrder(PetscSpace, PetscInt);
64: PETSC_EXTERN PetscErrorCode PetscSpaceGetOrder(PetscSpace, PetscInt *);
65: PETSC_EXTERN PetscErrorCode PetscSpaceEvaluate(PetscSpace, PetscInt, const PetscReal[], PetscReal[], PetscReal[], PetscReal[]);
66: PETSC_EXTERN PetscErrorCode PetscSpaceGetHeightSubspace(PetscSpace, PetscInt, PetscSpace *);
68: PETSC_EXTERN PetscErrorCode PetscSpacePolynomialSetNumVariables(PetscSpace, PetscInt);
69: PETSC_EXTERN PetscErrorCode PetscSpacePolynomialGetNumVariables(PetscSpace, PetscInt *);
70: PETSC_EXTERN PetscErrorCode PetscSpacePolynomialSetSymmetric(PetscSpace, PetscBool);
71: PETSC_EXTERN PetscErrorCode PetscSpacePolynomialGetSymmetric(PetscSpace, PetscBool *);
72: PETSC_EXTERN PetscErrorCode PetscSpacePolynomialSetTensor(PetscSpace, PetscBool);
73: PETSC_EXTERN PetscErrorCode PetscSpacePolynomialGetTensor(PetscSpace, PetscBool *);
75: PETSC_EXTERN PetscErrorCode PetscSpacePointGetPoints(PetscSpace, PetscQuadrature *);
76: PETSC_EXTERN PetscErrorCode PetscSpacePointSetPoints(PetscSpace, PetscQuadrature);
78: PETSC_EXTERN PetscClassId PETSCDUALSPACE_CLASSID;
80: /*J
81: PetscDualSpaceType - String with the name of a PETSc dual space
83: Level: beginner
85: .seealso: PetscDualSpaceSetType(), PetscDualSpace
86: J*/
87: typedef const char *PetscDualSpaceType;
88: #define PETSCDUALSPACELAGRANGE "lagrange"
89: #define PETSCDUALSPACESIMPLE "simple"
91: PETSC_EXTERN PetscFunctionList PetscDualSpaceList;
92: PETSC_EXTERN PetscErrorCode PetscDualSpaceCreate(MPI_Comm, PetscDualSpace *);
93: PETSC_EXTERN PetscErrorCode PetscDualSpaceDestroy(PetscDualSpace *);
94: PETSC_EXTERN PetscErrorCode PetscDualSpaceDuplicate(PetscDualSpace, PetscDualSpace *);
95: PETSC_EXTERN PetscErrorCode PetscDualSpaceSetType(PetscDualSpace, PetscDualSpaceType);
96: PETSC_EXTERN PetscErrorCode PetscDualSpaceGetType(PetscDualSpace, PetscDualSpaceType *);
97: PETSC_EXTERN PetscErrorCode PetscDualSpaceGetNumDof(PetscDualSpace, const PetscInt **);
98: PETSC_EXTERN PetscErrorCode PetscDualSpaceSetUp(PetscDualSpace);
99: PETSC_EXTERN PetscErrorCode PetscDualSpaceSetFromOptions(PetscDualSpace);
100: PETSC_STATIC_INLINE PetscErrorCode PetscDualSpaceViewFromOptions(PetscDualSpace A,PetscObject B,const char name[]) {return PetscObjectViewFromOptions((PetscObject)A,B,name);}
102: PETSC_EXTERN PetscErrorCode PetscDualSpaceView(PetscDualSpace,PetscViewer);
103: PETSC_EXTERN PetscErrorCode PetscDualSpaceRegister(const char [], PetscErrorCode (*)(PetscDualSpace));
104: PETSC_EXTERN PetscErrorCode PetscDualSpaceRegisterDestroy(void);
106: PETSC_EXTERN PetscErrorCode PetscDualSpaceGetDimension(PetscDualSpace, PetscInt *);
107: PETSC_EXTERN PetscErrorCode PetscDualSpaceSetNumComponents(PetscDualSpace, PetscInt);
108: PETSC_EXTERN PetscErrorCode PetscDualSpaceGetNumComponents(PetscDualSpace, PetscInt *);
109: PETSC_EXTERN PetscErrorCode PetscDualSpaceSetOrder(PetscDualSpace, PetscInt);
110: PETSC_EXTERN PetscErrorCode PetscDualSpaceGetOrder(PetscDualSpace, PetscInt *);
111: PETSC_EXTERN PetscErrorCode PetscDualSpaceSetDM(PetscDualSpace, DM);
112: PETSC_EXTERN PetscErrorCode PetscDualSpaceGetDM(PetscDualSpace, DM *);
113: PETSC_EXTERN PetscErrorCode PetscDualSpaceGetFunctional(PetscDualSpace, PetscInt, PetscQuadrature *);
114: PETSC_EXTERN PetscErrorCode PetscDualSpaceCreateReferenceCell(PetscDualSpace, PetscInt, PetscBool, DM *);
115: PETSC_EXTERN PetscErrorCode PetscDualSpaceGetSymmetries(PetscDualSpace, const PetscInt ****, const PetscScalar ****);
117: PETSC_EXTERN PetscErrorCode PetscDualSpaceApply(PetscDualSpace, PetscInt, PetscReal, PetscFECellGeom *, PetscInt, PetscErrorCode (*)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void *, PetscScalar *);
118: PETSC_EXTERN PetscErrorCode PetscDualSpaceApplyDefault(PetscDualSpace, PetscInt, PetscReal, PetscFECellGeom *, PetscInt, PetscErrorCode (*)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void *, PetscScalar *);
120: PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeGetContinuity(PetscDualSpace, PetscBool *);
121: PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeSetContinuity(PetscDualSpace, PetscBool);
122: PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeGetTensor(PetscDualSpace, PetscBool *);
123: PETSC_EXTERN PetscErrorCode PetscDualSpaceLagrangeSetTensor(PetscDualSpace, PetscBool);
125: PETSC_EXTERN PetscErrorCode PetscDualSpaceGetHeightSubspace(PetscDualSpace,PetscInt,PetscDualSpace *);
126: PETSC_EXTERN PetscErrorCode PetscDualSpaceSimpleSetDimension(PetscDualSpace, PetscInt);
127: PETSC_EXTERN PetscErrorCode PetscDualSpaceSimpleSetFunctional(PetscDualSpace, PetscInt, PetscQuadrature);
129: PETSC_EXTERN PetscClassId PETSCFE_CLASSID;
131: /*J
132: PetscFEType - String with the name of a PETSc finite element space
134: Level: beginner
136: Note: Currently, the classes are concerned with the implementation of element integration
138: .seealso: PetscFESetType(), PetscFE
139: J*/
140: typedef const char *PetscFEType;
141: #define PETSCFEBASIC "basic"
142: #define PETSCFENONAFFINE "nonaffine"
143: #define PETSCFEOPENCL "opencl"
144: #define PETSCFECOMPOSITE "composite"
146: PETSC_EXTERN PetscFunctionList PetscFEList;
147: PETSC_EXTERN PetscErrorCode PetscFECreate(MPI_Comm, PetscFE *);
148: PETSC_EXTERN PetscErrorCode PetscFEDestroy(PetscFE *);
149: PETSC_EXTERN PetscErrorCode PetscFESetType(PetscFE, PetscFEType);
150: PETSC_EXTERN PetscErrorCode PetscFEGetType(PetscFE, PetscFEType *);
151: PETSC_EXTERN PetscErrorCode PetscFESetUp(PetscFE);
152: PETSC_EXTERN PetscErrorCode PetscFESetFromOptions(PetscFE);
153: PETSC_STATIC_INLINE PetscErrorCode PetscFEViewFromOptions(PetscFE A,PetscObject B,const char name[]) {return PetscObjectViewFromOptions((PetscObject)A,B,name);}
155: PETSC_EXTERN PetscErrorCode PetscFEView(PetscFE,PetscViewer);
156: PETSC_EXTERN PetscErrorCode PetscFERegister(const char [], PetscErrorCode (*)(PetscFE));
157: PETSC_EXTERN PetscErrorCode PetscFERegisterDestroy(void);
158: PETSC_EXTERN PetscErrorCode PetscFECreateDefault(DM, PetscInt, PetscInt, PetscBool, const char [], PetscInt, PetscFE *);
160: PETSC_EXTERN PetscErrorCode PetscFEGetDimension(PetscFE, PetscInt *);
161: PETSC_EXTERN PetscErrorCode PetscFEGetSpatialDimension(PetscFE, PetscInt *);
162: PETSC_EXTERN PetscErrorCode PetscFESetNumComponents(PetscFE, PetscInt);
163: PETSC_EXTERN PetscErrorCode PetscFEGetNumComponents(PetscFE, PetscInt *);
164: PETSC_EXTERN PetscErrorCode PetscFEGetTileSizes(PetscFE, PetscInt *, PetscInt *, PetscInt *, PetscInt *);
165: PETSC_EXTERN PetscErrorCode PetscFESetTileSizes(PetscFE, PetscInt, PetscInt, PetscInt, PetscInt);
166: PETSC_EXTERN PetscErrorCode PetscFESetBasisSpace(PetscFE, PetscSpace);
167: PETSC_EXTERN PetscErrorCode PetscFEGetBasisSpace(PetscFE, PetscSpace *);
168: PETSC_EXTERN PetscErrorCode PetscFESetDualSpace(PetscFE, PetscDualSpace);
169: PETSC_EXTERN PetscErrorCode PetscFEGetDualSpace(PetscFE, PetscDualSpace *);
170: PETSC_EXTERN PetscErrorCode PetscFESetQuadrature(PetscFE, PetscQuadrature);
171: PETSC_EXTERN PetscErrorCode PetscFEGetQuadrature(PetscFE, PetscQuadrature *);
172: PETSC_EXTERN PetscErrorCode PetscFESetFaceQuadrature(PetscFE, PetscQuadrature);
173: PETSC_EXTERN PetscErrorCode PetscFEGetFaceQuadrature(PetscFE, PetscQuadrature *);
174: PETSC_EXTERN PetscErrorCode PetscFEGetNumDof(PetscFE, const PetscInt **);
175: PETSC_EXTERN PetscErrorCode PetscFEGetDefaultTabulation(PetscFE, PetscReal **, PetscReal **, PetscReal **);
176: PETSC_EXTERN PetscErrorCode PetscFEGetFaceTabulation(PetscFE, PetscReal **, PetscReal **, PetscReal **);
177: PETSC_EXTERN PetscErrorCode PetscFEGetFaceCentroidTabulation(PetscFE, PetscReal **);
178: PETSC_EXTERN PetscErrorCode PetscFEGetTabulation(PetscFE, PetscInt, const PetscReal[], PetscReal **, PetscReal **, PetscReal **);
179: PETSC_EXTERN PetscErrorCode PetscFERestoreTabulation(PetscFE, PetscInt, const PetscReal[], PetscReal **, PetscReal **, PetscReal **);
180: PETSC_EXTERN PetscErrorCode PetscFERefine(PetscFE, PetscFE *);
181: PETSC_EXTERN PetscErrorCode PetscFEGetHeightSubspace(PetscFE, PetscInt, PetscFE *);
183: PETSC_EXTERN PetscErrorCode PetscFEIntegrate(PetscFE, PetscDS, PetscInt, PetscInt, PetscFECellGeom *, const PetscScalar[], PetscDS, const PetscScalar[], PetscScalar[]);
184: PETSC_EXTERN PetscErrorCode PetscFEIntegrateResidual(PetscFE, PetscDS, PetscInt, PetscInt, PetscFECellGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscScalar[]);
185: PETSC_EXTERN PetscErrorCode PetscFEIntegrateBdResidual(PetscFE, PetscDS, PetscInt, PetscInt, PetscFEFaceGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscScalar[]);
186: PETSC_EXTERN PetscErrorCode PetscFEIntegrateJacobian(PetscFE, PetscDS, PetscFEJacobianType, PetscInt, PetscInt, PetscInt, PetscFECellGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscReal, PetscScalar[]);
187: PETSC_EXTERN PetscErrorCode PetscFEIntegrateBdJacobian(PetscFE, PetscDS, PetscInt, PetscInt, PetscInt, PetscFEFaceGeom *, const PetscScalar[], const PetscScalar[], PetscDS, const PetscScalar[], PetscReal, PetscReal, PetscScalar[]);
189: PETSC_EXTERN PetscErrorCode PetscFECompositeGetMapping(PetscFE, PetscInt *, const PetscReal *[], const PetscReal *[], const PetscReal *[]);
191: PETSC_EXTERN PetscErrorCode PetscFEOpenCLSetRealType(PetscFE, PetscDataType);
192: PETSC_EXTERN PetscErrorCode PetscFEOpenCLGetRealType(PetscFE, PetscDataType *);
194: #endif