Basis
libCEED internally uses row-major (C convention) storage of matrices, while Julia uses column-major (Fortran convention) storage.
LibCEED.jl will typically handle the conversion between these formats by transposing or permuting the dimensions of the input and output matrices and tensors.
LibCEED.Basis
— TypeBasis
Wraps a CeedBasis
object, representing a finite element basis. A Basis
object can be created using one of:
Missing docstring for BasisCollocated
. Check Documenter's build log for details.
LibCEED.create_tensor_h1_lagrange_basis
— Functioncreate_tensor_h1_lagrange_basis(ceed, dim, ncomp, p, q, qmode)
Create a tensor-product Lagrange basis.
Arguments:
ceed
: ACeed
object where theBasis
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).
LibCEED.create_tensor_h1_basis
— Functioncreate_tensor_h1_basis(c::Ceed, dim, ncomp, p, q, interp1d, grad1d, qref1d, qweight1d)
Create a tensor-product basis for $H^1$ discretizations.
Arguments:
ceed
: ACeed
object where theBasis
will be created.dim
: Topological dimension.ncomp
: Number of field components (1 for scalar fields).p
: Number of nodes in one dimension.q
: Number of quadrature points in one dimensioninterp1d
: Matrix of size(q, p)
expressing the values of nodal basis functions at quadrature points.grad1d
: Matrix of size(p, q)
expressing derivatives of nodal basis functions at quadrature points.qref1d
: Array of lengthq
holding the locations of quadrature points on the 1D reference element $[-1, 1]$.qweight1d
: Array of lengthq
holding the quadrature weights on the reference element.
LibCEED.create_h1_basis
— Functioncreate_h1_basis(c::Ceed, topo::Topology, ncomp, nnodes, nqpts, interp, grad, qref, qweight)
Create a non tensor-product basis for $H^1$ discretizations
Arguments:
ceed
: ACeed
object where theBasis
will be created.topo
:Topology
of element, e.g. hypercube, simplex, etc.ncomp
: Number of field components (1 for scalar fields).nnodes
: Total number of nodes.nqpts
: Total number of quadrature points.interp
: Matrix of size(nqpts, nnodes)
expressing the values of nodal basis functions at quadrature points.grad
: Array of size(dim, nqpts, nnodes)
expressing derivatives of nodal basis functions at quadrature points.qref
: Matrix of size(dim, nqpts)
holding the locations of quadrature points on the reference element $[-1, 1]$.qweight
: Array of lengthnqpts
holding the quadrature weights on the reference element.
LibCEED.create_hdiv_basis
— Functioncreate_hdiv_basis(c::Ceed, topo::Topology, ncomp, nnodes, nqpts, interp, div, qref, qweight)
Create a non tensor-product basis for H(div) discretizations
Arguments:
ceed
: ACeed
object where theBasis
will be created.topo
:Topology
of element, e.g. hypercube, simplex, etc.ncomp
: Number of field components (1 for scalar fields).nnodes
: Total number of nodes.nqpts
: Total number of quadrature points.interp
: Matrix of size(dim, nqpts, nnodes)
expressing the values of basis functions at quadrature points.div
: Array of size(nqpts, nnodes)
expressing divergence of basis functions at quadrature points.qref
: Matrix of size(dim, nqpts)
holding the locations of quadrature points on the reference element $[-1, 1]$.qweight
: Array of lengthnqpts
holding the quadrature weights on the reference element.
LibCEED.create_hcurl_basis
— Functioncreate_hcurl_basis(c::Ceed, topo::Topology, ncomp, nnodes, nqpts, interp, curl, qref, qweight)
Create a non tensor-product basis for H(curl) discretizations
Arguments:
ceed
: ACeed
object where theBasis
will be created.topo
:Topology
of element, e.g. hypercube, simplex, etc.ncomp
: Number of field components (1 for scalar fields).nnodes
: Total number of nodes.nqpts
: Total number of quadrature points.interp
: Matrix of size(dim, nqpts, nnodes)
expressing the values of basis functions at quadrature points.curl
: Matrix of size(curlcomp, nqpts, nnodes)
,curlcomp = 1 if dim < 3 else dim
) matrix expressing curl of basis functions at quadrature points.qref
: Matrix of size(dim, nqpts)
holding the locations of quadrature points on the reference element $[-1, 1]$.qweight
: Array of lengthnqpts
holding the quadrature weights on the reference element.
LibCEED.apply!
— Methodapply!(b::Basis, nelem, tmode::TransposeMode, emode::EvalMode, u::AbstractCeedVector, v::AbstractCeedVector)
Apply basis evaluation from nodes to quadrature points or vice versa, storing the result in the CeedVector
v
.
nelem
specifies the number of elements to apply the basis evaluation to; the backend will specify the ordering in CeedElemRestrictionCreateBlocked()
Set tmode
to CEED_NOTRANSPOSE
to evaluate from nodes to quadrature or to CEED_TRANSPOSE
to apply the transpose, mapping from quadrature points to nodes.
Set the EvalMode
emode
to:
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.
LibCEED.apply
— Methodapply(b::Basis, u::AbstractVector; nelem=1, tmode=NOTRANSPOSE, emode=EVAL_INTERP)
Performs the same function as the above-defined apply!
, but automatically convert from Julia arrays to CeedVector
for convenience.
The result will be returned in a newly allocated array of the correct size.
LibCEED.getdimension
— Functiongetdimension(b::Basis)
Return the spatial dimension of the given Basis
.
getdimension(t::Topology)
Return the spatial dimension of the given Topology
.
LibCEED.gettopology
— FunctionLibCEED.getnumcomponents
— Functiongetnumcomponents(b::Basis)
Return the number of components of the given Basis
.
getnumcomponents(r::ElemRestriction)
Get the number of components in the elements of an ElemRestriction
.
LibCEED.getnumnodes
— Functiongetnumnodes(b::Basis)
Return the number of nodes of the given Basis
.
LibCEED.getnumnodes1d
— Functiongetnumnodes1d(b::Basis)
Return the number of 1D nodes of the given (tensor-product) [`Basis`](@ref).
LibCEED.getnumqpts
— Functiongetnumqpts(b::Basis)
Return the number of quadrature points of the given Basis
.
LibCEED.getnumqpts1d
— Functiongetnumqpts1d(b::Basis)
Return the number of 1D quadrature points of the given (tensor-product) Basis
.
LibCEED.getqref
— Functiongetqref(b::Basis)
Get the reference coordinates of quadrature points (in dim
dimensions) of the given Basis
.
LibCEED.getqweights
— Functiongetqref(b::Basis)
Get the quadrature weights of quadrature points (in dim
dimensions) of the given Basis
.
LibCEED.getinterp
— Functiongetinterp(b::Basis)
Get the interpolation matrix of the given Basis
. Returns a matrix of size (getnumqpts(b), getnumnodes(b))
for a given $H^1$ basis or (getdimension(b), getnumqpts(b), getnumnodes(b))
for a given vector $H(div)$ or $H(curl)$ basis.
LibCEED.getinterp1d
— Functiongetinterp1d(b::Basis)
Get the 1D interpolation matrix of the given Basis
. b
must be a tensor-product basis, otherwise this function will fail. Returns a matrix of size (getnumqpts1d(b), getnumnodes1d(b))
.
LibCEED.getgrad
— Functiongetgrad(b::Basis)
Get the gradient matrix of the given Basis
. Returns a tensor of size (getdimension(b), getnumqpts(b), getnumnodes(b))
.
LibCEED.getgrad1d
— Functiongetgrad1d(b::Basis)
Get the 1D derivative matrix of the given Basis
. Returns a matrix of size (getnumqpts(b), getnumnodes(b))
.
Missing docstring for getdiv
. Check Documenter's build log for details.
Missing docstring for getcurl
. Check Documenter's build log for details.