CeedBasis¶
A CeedBasis defines the discrete finite element basis and associated quadrature rule.
Discrete element bases and quadrature¶
-
typedef struct CeedBasis_private *
CeedBasis
¶ Handle for object describing discrete finite element evaluations.
-
const CeedBasis
CEED_BASIS_COLLOCATED
= &ceed_basis_collocated¶ Indicate that the quadrature points are collocated with the nodes.
-
int
CeedBasisCreateTensorH1
(Ceed ceed, CeedInt dim, CeedInt ncomp, CeedInt P1d, CeedInt Q1d, const CeedScalar *interp1d, const CeedScalar *grad1d, const CeedScalar *qref1d, const CeedScalar *qweight1d, CeedBasis *basis)¶ Create a tensor-product basis for H^1 discretizations.
User Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
ceed
: A Ceed object where the CeedBasis will be createddim
: Topological dimensionncomp
: Number of field components (1 for scalar fields)P1d
: Number of nodes in one dimensionQ1d
: Number of quadrature points in one dimensioninterp1d
: Row-major (Q1d * P1d) matrix expressing the values of nodal basis functions at quadrature pointsgrad1d
: Row-major (Q1d * P1d) matrix expressing derivatives of nodal basis functions at quadrature pointsqref1d
: Array of length Q1d holding the locations of quadrature points on the 1D reference element [-1, 1]qweight1d
: Array of length Q1d holding the quadrature weights on the reference element[out] basis
: Address of the variable where the newly created CeedBasis will be stored.
-
int
CeedBasisCreateTensorH1Lagrange
(Ceed ceed, CeedInt dim, CeedInt ncomp, CeedInt P, CeedInt Q, CeedQuadMode qmode, CeedBasis *basis)¶ Create a tensor-product Lagrange basis.
User Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
ceed
: A Ceed object where the CeedBasis will be createddim
: Topological dimension of elementncomp
: Number of field components (1 for scalar fields)P
: Number of Gauss-Lobatto nodes in one dimension. The polynomial degree of the resulting Q_k element is k=P-1.Q
: Number of quadrature points in one dimension.qmode
: Distribution of the Q quadrature points (affects order of accuracy for the quadrature)[out] basis
: Address of the variable where the newly created CeedBasis will be stored.
-
int
CeedBasisCreateH1
(Ceed ceed, CeedElemTopology topo, CeedInt ncomp, CeedInt nnodes, CeedInt nqpts, const CeedScalar *interp, const CeedScalar *grad, const CeedScalar *qref, const CeedScalar *qweight, CeedBasis *basis)¶ Create a non tensor-product basis for H^1 discretizations.
User Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
ceed
: A Ceed object where the CeedBasis will be createdtopo
: Topology of element, e.g. hypercube, simplex, ectncomp
: Number of field components (1 for scalar fields)nnodes
: Total number of nodesnqpts
: Total number of quadrature pointsinterp
: Row-major (nqpts * nnodes) matrix expressing the values of nodal basis functions at quadrature pointsgrad
: Row-major (nqpts * dim * nnodes) matrix expressing derivatives of nodal basis functions at quadrature pointsqref
: Array of length nqpts holding the locations of quadrature points on the reference element [-1, 1]qweight
: Array of length nqpts holding the quadrature weights on the reference element[out] basis
: Address of the variable where the newly created CeedBasis will be stored.
-
int
CeedBasisView
(CeedBasis basis, FILE *stream)¶ View a CeedBasis.
User Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis to viewstream
: Stream to view to, e.g., stdout
-
int
CeedBasisApply
(CeedBasis basis, CeedInt nelem, CeedTransposeMode tmode, CeedEvalMode emode, CeedVector u, CeedVector v)¶ Apply basis evaluation from nodes to quadrature points or vice versa.
User Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis to evaluatenelem
: The number of elements to apply the basis evaluation to; the backend will specify the ordering in CeedElemRestrictionCreateBlocked()tmode
: CEED_NOTRANSPOSE to evaluate from nodes to quadrature points, CEED_TRANSPOSE to apply the transpose, mapping from quadrature points to nodesemode
: CEED_EVAL_NONE to use values directly, CEED_EVAL_INTERP to use interpolated values, CEED_EVAL_GRAD to use gradients, CEED_EVAL_WEIGHT to use quadrature weights.[in] u
: Input CeedVector[out] v
: Output CeedVector
-
int
CeedBasisGetDimension
(CeedBasis basis, CeedInt *dim)¶ Get dimension for given CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] dim
: Variable to store dimension of basis
-
int
CeedBasisGetTopology
(CeedBasis basis, CeedElemTopology *topo)¶ Get topology for given CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] topo
: Variable to store topology of basis
-
int
CeedBasisGetNumComponents
(CeedBasis basis, CeedInt *numcomp)¶ Get number of components for given CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] numcomp
: Variable to store number of components of basis
-
int
CeedBasisGetNumNodes
(CeedBasis basis, CeedInt *P)¶ Get total number of nodes (in dim dimensions) of a CeedBasis.
Utility Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] P
: Variable to store number of nodes
-
int
CeedBasisGetNumNodes1D
(CeedBasis basis, CeedInt *P1d)¶ Get total number of nodes (in 1 dimension) of a CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] P1d
: Variable to store number of nodes
-
int
CeedBasisGetNumQuadraturePoints
(CeedBasis basis, CeedInt *Q)¶ Get total number of quadrature points (in dim dimensions) of a CeedBasis.
Utility Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] Q
: Variable to store number of quadrature points
-
int
CeedBasisGetNumQuadraturePoints1D
(CeedBasis basis, CeedInt *Q1d)¶ Get total number of quadrature points (in 1 dimension) of a CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] Q1d
: Variable to store number of quadrature points
-
int
CeedBasisGetQRef
(CeedBasis basis, const CeedScalar **qref)¶ Get reference coordinates of quadrature points (in dim dimensions) of a CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] qref
: Variable to store reference coordinates of quadrature points
-
int
CeedBasisGetQWeights
(CeedBasis basis, const CeedScalar **qweight)¶ Get quadrature weights of quadrature points (in dim dimensions) of a CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] qweight
: Variable to store quadrature weights
-
int
CeedBasisGetInterp
(CeedBasis basis, const CeedScalar **interp)¶ Get interpolation matrix of a CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] interp
: Variable to store interpolation matrix
-
int
CeedBasisGetInterp1D
(CeedBasis basis, const CeedScalar **interp1d)¶ Get 1D interpolation matrix of a tensor product CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] interp1d
: Variable to store interpolation matrix
-
int
CeedBasisGetGrad
(CeedBasis basis, const CeedScalar **grad)¶ Get gradient matrix of a CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] grad
: Variable to store gradient matrix
-
int
CeedBasisGetGrad1D
(CeedBasis basis, const CeedScalar **grad1d)¶ Get 1D gradient matrix of a tensor product CeedBasis.
Backend Developer Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis[out] grad1d
: Variable to store gradient matrix
-
int
CeedBasisDestroy
(CeedBasis *basis)¶ Destroy a CeedBasis.
User Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
basis
: CeedBasis to destroy
-
int
CeedGaussQuadrature
(CeedInt Q, CeedScalar *qref1d, CeedScalar *qweight1d)¶ Construct a Gauss-Legendre quadrature.
Utility Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
Q
: Number of quadrature points (integrates polynomials of degree 2*Q-1 exactly)[out] qref1d
: Array of length Q to hold the abscissa on [-1, 1][out] qweight1d
: Array of length Q to hold the weights
-
int
CeedLobattoQuadrature
(CeedInt Q, CeedScalar *qref1d, CeedScalar *qweight1d)¶ Construct a Gauss-Legendre-Lobatto quadrature.
Utility Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
Q
: Number of quadrature points (integrates polynomials of degree 2*Q-3 exactly)[out] qref1d
: Array of length Q to hold the abscissa on [-1, 1][out] qweight1d
: Array of length Q to hold the weights
-
int
CeedQRFactorization
(Ceed ceed, CeedScalar *mat, CeedScalar *tau, CeedInt m, CeedInt n)¶ Return QR Factorization of a matrix.
Utility Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
ceed
: A Ceed context for error handling[inout] mat
: Row-major matrix to be factorized in place[inout] tau
: Vector of length m of scaling factorsm
: Number of rowsn
: Number of columns
-
int
CeedSymmetricSchurDecomposition
(Ceed ceed, CeedScalar *mat, CeedScalar *lambda, CeedInt n)¶ Return symmetric Schur decomposition of the symmetric matrix mat via symmetric QR factorization.
Utility Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
ceed
: A Ceed context for error handling[inout] mat
: Row-major matrix to be factorized in place[out] lambda
: Vector of length n of eigenvaluesn
: Number of rows/columns
-
int
CeedSimultaneousDiagonalization
(Ceed ceed, CeedScalar *matA, CeedScalar *matB, CeedScalar *x, CeedScalar *lambda, CeedInt n)¶ Return Simultaneous Diagonalization of two matrices.
This solves the generalized eigenvalue problem A x = lambda B x, where A and B are symmetric and B is positive definite. We generate the matrix X and vector Lambda such that X^T A X = Lambda and X^T B X = I. This is equivalent to the LAPACK routine ‘sygv’ with TYPE = 1.
Utility Functions
- Return
An error code: 0 - success, otherwise - failure
- Parameters
ceed
: A Ceed context for error handling[in] matA
: Row-major matrix to be factorized with eigenvalues[in] matB
: Row-major matrix to be factorized to identity[out] x
: Row-major orthogonal matrix[out] lambda
: Vector of length n of generalized eigenvaluesn
: Number of rows/columns
Typedefs and Enumerations¶
-
enum
CeedTransposeMode
¶ Denotes whether a linear transformation or its transpose should be applied.
Values:
-
enumerator
CEED_NOTRANSPOSE
¶ Apply the linear transformation.
-
enumerator
CEED_TRANSPOSE
¶ Apply the transpose.
-
enumerator
-
enum
CeedEvalMode
¶ Basis evaluation mode.
Modes can be bitwise ORed when passing to most functions.
Values:
-
enumerator
CEED_EVAL_NONE
¶ Perform no evaluation (either because there is no data or it is already at quadrature points)
-
enumerator
CEED_EVAL_INTERP
¶ Interpolate from nodes to quadrature points.
-
enumerator
CEED_EVAL_GRAD
¶ Evaluate gradients at quadrature points from input in a nodal basis.
-
enumerator
CEED_EVAL_DIV
¶ Evaluate divergence at quadrature points from input in a nodal basis.
-
enumerator
CEED_EVAL_CURL
¶ Evaluate curl at quadrature points from input in a nodal basis.
-
enumerator
CEED_EVAL_WEIGHT
¶ Using no input, evaluate quadrature weights on the reference element.
-
enumerator
-
enum
CeedQuadMode
¶ Type of quadrature; also used for location of nodes.
Values:
-
enumerator
CEED_GAUSS
¶ Gauss-Legendre quadrature.
-
enumerator
CEED_GAUSS_LOBATTO
¶ Gauss-Legendre-Lobatto quadrature.
-
enumerator
-
enum
CeedElemTopology
¶ Type of basis shape to create non-tensor H1 element basis.
Dimension can be extracted with bitwise AND (CeedElemTopology & 2**(dim + 2)) == TRUE
Values:
-
enumerator
CEED_LINE
¶ Line.
-
enumerator
CEED_TRIANGLE
¶ Triangle - 2D shape.
-
enumerator
CEED_QUAD
¶ Quadralateral - 2D shape.
-
enumerator
CEED_TET
¶ Tetrahedron - 3D shape.
-
enumerator
CEED_PYRAMID
¶ Pyramid - 3D shape.
-
enumerator
CEED_PRISM
¶ Prism - 3D shape.
-
enumerator
CEED_HEX
¶ Hexehedron - 3D shape.
-
enumerator