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 created

  • dim: Topological dimension

  • ncomp: Number of field components (1 for scalar fields)

  • P1d: Number of nodes in one dimension

  • Q1d: Number of quadrature points in one dimension

  • interp1d: Row-major (Q1d * P1d) matrix expressing the values of nodal basis functions at quadrature points

  • grad1d: Row-major (Q1d * P1d) matrix expressing derivatives of nodal basis functions at quadrature points

  • qref1d: 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 created

  • dim: Topological dimension of element

  • ncomp: 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 created

  • topo: Topology of element, e.g. hypercube, simplex, ect

  • ncomp: Number of field components (1 for scalar fields)

  • nnodes: Total number of nodes

  • nqpts: Total number of quadrature points

  • interp: Row-major (nqpts * nnodes) matrix expressing the values of nodal basis functions at quadrature points

  • grad: Row-major (nqpts * dim * nnodes) matrix expressing derivatives of nodal basis functions at quadrature points

  • qref: 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 view

  • stream: 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

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 factors

  • m: Number of rows

  • n: 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 eigenvalues

  • n: 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 eigenvalues

  • n: 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.

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.

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.

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.