Actual source code: petscdmceed.h
1: #pragma once
3: #include <petscdm.h>
5: /* MANSEC = DM */
7: #if defined(PETSC_HAVE_LIBCEED)
8: #include <ceed.h>
10: #if defined(PETSC_CLANG_STATIC_ANALYZER)
11: void PetscCallCEED(CeedErrorType);
12: #else
13: #define PetscCallCEED(...) \
14: do { \
15: CeedErrorType ierr_ceed_ = __VA_ARGS__; \
16: PetscCheck(ierr_ceed_ == CEED_ERROR_SUCCESS, PETSC_COMM_SELF, PETSC_ERR_LIB, "libCEED error: %s", CeedErrorTypes[ierr_ceed_]); \
17: } while (0)
18: #endif /* PETSC_CLANG_STATIC_ANALYZER */
19: #define CHKERRQ_CEED(...) PetscCallCEED(__VA_ARGS__)
21: PETSC_EXTERN PetscErrorCode DMGetCeed(DM, Ceed *);
23: PETSC_EXTERN PetscErrorCode VecGetCeedVector(Vec, Ceed, CeedVector *);
24: PETSC_EXTERN PetscErrorCode VecGetCeedVectorRead(Vec, Ceed, CeedVector *);
25: PETSC_EXTERN PetscErrorCode VecRestoreCeedVector(Vec, CeedVector *);
26: PETSC_EXTERN PetscErrorCode VecRestoreCeedVectorRead(Vec, CeedVector *);
27: PETSC_INTERN PetscErrorCode DMCeedCreate_Internal(DM, IS, PetscBool, CeedQFunctionUser, const char *, DMCeed *);
28: PETSC_EXTERN PetscErrorCode DMCeedCreate(DM, PetscBool, CeedQFunctionUser, const char *);
29: PETSC_EXTERN PetscErrorCode DMCeedCreateFVM(DM, PetscBool, CeedQFunctionUser, const char *, CeedQFunctionContext);
31: struct _PETSc_DMCEED {
32: CeedBasis basis; // Basis for element function space
33: CeedElemRestriction er; // Map from PETSc local vector to FE element vectors
34: CeedElemRestriction erL; // Map from PETSc local vector to FV left cell vectors
35: CeedElemRestriction erR; // Map from PETSc local vector to FV right cell vectors
36: CeedQFunctionUser func; // Plex Function for this operator
37: char *funcSource; // Plex Function source as text
38: CeedQFunction qf; // QFunction expressing the operator action
39: CeedOperator op; // Operator action for this object
40: DMCeed geom; // Operator computing geometric data at quadrature points
41: CeedElemRestriction erq; // Map from PETSc local vector to quadrature points
42: CeedVector qd; // Geometric data at quadrature points used in calculating the qfunction
43: DMCeed info; // Mesh information at quadrature points
44: CeedElemRestriction eri; // Map from PETSc local vector to quadrature points
45: CeedVector qi; // Mesh information at quadrature points
46: };
48: #else
50: struct _PETSc_DMCEED {
51: PetscInt dummy;
52: };
54: #endif
56: PETSC_EXTERN PetscErrorCode DMCeedComputeGeometry(DM, DMCeed);
57: PETSC_EXTERN PetscErrorCode DMCeedComputeInfo(DM, DMCeed);
58: PETSC_EXTERN PetscErrorCode DMCeedDestroy(DMCeed *);