![]() |
MFEM v4.8.0
Finite element discretization library
|
#include <fe_l2.hpp>
Public Member Functions | |
L2_FuentesPyramidElement (const int p, const int btype=BasisType::GaussLegendre) | |
Construct the L2_PyramidElement of order p and BasisType btype. | |
virtual void | CalcShape (const IntegrationPoint &ip, Vector &shape) const |
Evaluate the values of all shape functions of a scalar finite element in reference space at the given point ip. | |
virtual void | CalcDShape (const IntegrationPoint &ip, DenseMatrix &dshape) const |
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the given point ip. | |
![]() | |
NodalFiniteElement (int D, Geometry::Type G, int Do, int O, int F=FunctionSpace::Pk) | |
Construct NodalFiniteElement with given. | |
const DofToQuad & | GetDofToQuad (const IntegrationRule &ir, DofToQuad::Mode mode) const override |
Return a DofToQuad structure corresponding to the given IntegrationRule using the given DofToQuad::Mode. | |
void | GetLocalInterpolation (ElementTransformation &Trans, DenseMatrix &I) const override |
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base geometry under the given transformation. | |
void | GetLocalRestriction (ElementTransformation &Trans, DenseMatrix &R) const override |
Return a local restriction matrix R (Dof x Dof) mapping fine dofs to coarse dofs. | |
void | GetTransferMatrix (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const override |
Return interpolation matrix, I, which maps dofs from a coarse element, fe, to the fine dofs on this finite element. | |
void | Project (Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const override |
Given a coefficient and a transformation, compute its projection (approximation) in the local finite dimensional space in terms of the degrees of freedom. | |
void | Project (VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override |
Given a vector coefficient and a transformation, compute its projection (approximation) in the local finite dimensional space in terms of the degrees of freedom. (VectorFiniteElements) | |
void | ProjectMatrixCoefficient (MatrixCoefficient &mc, ElementTransformation &T, Vector &dofs) const override |
Given a matrix coefficient and a transformation, compute an approximation ("projection") in the local finite dimensional space in terms of the degrees of freedom. For VectorFiniteElements, the rows of the coefficient are projected in the vector space. | |
void | Project (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const override |
Compute the embedding/projection matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the projection depends on it. | |
void | ProjectGrad (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const override |
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the matrix depends on it. | |
void | ProjectDiv (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &div) const override |
Compute the discrete divergence matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the matrix depends on it. | |
const Array< int > & | GetLexicographicOrdering () const |
Get an Array<int> that maps lexicographically ordered indices to the indices of the respective nodes/dofs/basis functions. | |
void | ReorderLexToNative (int ncomp, Vector &dofs) const |
![]() | |
ScalarFiniteElement (int D, Geometry::Type G, int Do, int O, int F=FunctionSpace::Pk) | |
Construct ScalarFiniteElement with given. | |
virtual void | SetMapType (int M) |
Set the FiniteElement::MapType of the element to either VALUE or INTEGRAL. Also sets the FiniteElement::DerivType to GRAD if the FiniteElement::MapType is VALUE. | |
void | NodalLocalInterpolation (ElementTransformation &Trans, DenseMatrix &I, const ScalarFiniteElement &fine_fe) const |
Get the matrix I that defines nodal interpolation between this element and the refined element fine_fe. | |
void | ScalarLocalInterpolation (ElementTransformation &Trans, DenseMatrix &I, const ScalarFiniteElement &fine_fe) const |
Get matrix I "Interpolation" defined through local L2-projection in the space defined by the fine_fe. | |
void | ScalarLocalL2Restriction (ElementTransformation &Trans, DenseMatrix &R, const ScalarFiniteElement &coarse_fe) const |
Get restriction matrix R defined through local L2-projection in the space defined by the coarse_fe. | |
![]() | |
FiniteElement (int D, Geometry::Type G, int Do, int O, int F=FunctionSpace::Pk) | |
Construct FiniteElement with given. | |
int | GetDim () const |
Returns the reference space dimension for the finite element. | |
int | GetRangeDim () const |
Returns the vector dimension for vector-valued finite elements, which is also the dimension of the interpolation operation. | |
int | GetCurlDim () const |
Returns the dimension of the curl for vector-valued finite elements. | |
Geometry::Type | GetGeomType () const |
Returns the Geometry::Type of the reference element. | |
int | GetDof () const |
Returns the number of degrees of freedom in the finite element. | |
int | GetOrder () const |
Returns the order of the finite element. In the case of anisotropic orders, returns the maximum order. | |
bool | HasAnisotropicOrders () const |
Returns true if the FiniteElement basis may be using different orders/degrees in different spatial directions. | |
const int * | GetAnisotropicOrders () const |
Returns an array containing the anisotropic orders/degrees. | |
int | Space () const |
Returns the type of FunctionSpace on the element. | |
int | GetRangeType () const |
Returns the FiniteElement::RangeType of the element, one of {SCALAR, VECTOR}. | |
int | GetDerivRangeType () const |
Returns the FiniteElement::RangeType of the element derivative, either SCALAR or VECTOR. | |
int | GetMapType () const |
Returns the FiniteElement::MapType of the element describing how reference functions are mapped to physical space, one of {VALUE, INTEGRAL H_DIV, H_CURL}. | |
int | GetDerivType () const |
Returns the FiniteElement::DerivType of the element describing the spatial derivative method implemented, one of {NONE, GRAD, DIV, CURL}. | |
int | GetDerivMapType () const |
Returns the FiniteElement::DerivType of the element describing how reference function derivatives are mapped to physical space, one of {VALUE, INTEGRAL, H_DIV, H_CURL}. | |
void | CalcPhysShape (ElementTransformation &Trans, Vector &shape) const |
Evaluate the values of all shape functions of a scalar finite element in physical space at the point described by Trans. | |
void | CalcPhysDShape (ElementTransformation &Trans, DenseMatrix &dshape) const |
Evaluate the gradients of all shape functions of a scalar finite element in physical space at the point described by Trans. | |
const IntegrationRule & | GetNodes () const |
Get a const reference to the nodes of the element. | |
virtual void | CalcHessian (const IntegrationPoint &ip, DenseMatrix &Hessian) const |
Evaluate the Hessians of all shape functions of a scalar finite element in reference space at the given point ip. | |
void | CalcPhysHessian (ElementTransformation &Trans, DenseMatrix &Hessian) const |
Evaluate the Hessian of all shape functions of a scalar finite element in physical space at the given point ip. | |
void | CalcPhysLaplacian (ElementTransformation &Trans, Vector &Laplacian) const |
Evaluate the Laplacian of all shape functions of a scalar finite element in physical space at the given point ip. | |
void | CalcPhysLinLaplacian (ElementTransformation &Trans, Vector &Laplacian) const |
Evaluate the Laplacian of all shape functions of a scalar finite element in physical space at the given point ip. | |
virtual void | CalcVShape (const IntegrationPoint &ip, DenseMatrix &shape) const |
Evaluate the values of all shape functions of a vector finite element in reference space at the given point ip. | |
virtual void | CalcVShape (ElementTransformation &Trans, DenseMatrix &shape) const |
Evaluate the values of all shape functions of a vector finite element in physical space at the point described by Trans. | |
void | CalcPhysVShape (ElementTransformation &Trans, DenseMatrix &shape) const |
Equivalent to the CalcVShape() method with the same arguments. | |
virtual void | CalcDivShape (const IntegrationPoint &ip, Vector &divshape) const |
Evaluate the divergence of all shape functions of a vector finite element in reference space at the given point ip. | |
void | CalcPhysDivShape (ElementTransformation &Trans, Vector &divshape) const |
Evaluate the divergence of all shape functions of a vector finite element in physical space at the point described by Trans. | |
virtual void | CalcCurlShape (const IntegrationPoint &ip, DenseMatrix &curl_shape) const |
Evaluate the curl of all shape functions of a vector finite element in reference space at the given point ip. | |
virtual void | CalcPhysCurlShape (ElementTransformation &Trans, DenseMatrix &curl_shape) const |
Evaluate the curl of all shape functions of a vector finite element in physical space at the point described by Trans. | |
virtual void | GetFaceDofs (int face, int **dofs, int *ndofs) const |
Get the dofs associated with the given face. *dofs is set to an internal array of the local dofc on the face, while *ndofs is set to the number of dofs on that face. | |
virtual void | ProjectFromNodes (Vector &vc, ElementTransformation &Trans, Vector &dofs) const |
Given a vector of values at the finite element nodes and a transformation, compute its projection (approximation) in the local finite dimensional space in terms of the degrees of freedom. Valid for VectorFiniteElements. | |
virtual void | ProjectDelta (int vertex, Vector &dofs) const |
Project a delta function centered on the given vertex in the local finite dimensional space represented by the dofs. | |
virtual void | ProjectCurl (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const |
Compute the discrete curl matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the matrix depends on it. | |
virtual void | GetFaceMap (const int face_id, Array< int > &face_map) const |
Return the mapping from lexicographic face DOFs to lexicographic element DOFs for the given local face face_id. | |
virtual const StatelessDofTransformation * | GetDofTransformation () const |
Return a DoF transformation object for this particular type of basis. | |
virtual | ~FiniteElement () |
Deconstruct the FiniteElement. | |
![]() | |
FuentesPyramid ()=default | |
void | phi_E (int p, Vector s, const DenseMatrix &grad_s, Vector &u, DenseMatrix &grad_u) const |
void | phi_Q (int p, Vector s, Vector t, DenseMatrix &u) const |
void | phi_Q (int p, Vector s, const DenseMatrix &grad_s, Vector t, const DenseMatrix &grad_t, DenseMatrix &u, DenseTensor &grad_u) const |
void | phi_T (int p, Vector nu, DenseMatrix &u) const |
void | phi_T (int p, Vector nu, const DenseMatrix &grad_nu, DenseMatrix &u, DenseTensor &grad_u) const |
void | E_E (int p, Vector s, Vector sds, DenseMatrix &u) const |
void | E_E (int p, Vector s, const DenseMatrix &grad_s, DenseMatrix &u, DenseMatrix &curl_u) const |
void | E_Q (int p, Vector s, Vector ds, Vector t, DenseTensor &u) const |
void | E_Q (int p, Vector s, const DenseMatrix &grad_s, Vector t, const DenseMatrix &grad_t, DenseTensor &u, DenseTensor &curl_u) const |
void | E_T (int p, Vector s, Vector sds, DenseTensor &u) const |
void | E_T (int p, Vector s, const DenseMatrix &grad_s, DenseTensor &u, DenseTensor &curl_u) const |
void | V_Q (int p, Vector s, Vector ds, Vector t, Vector dt, DenseTensor &u) const |
void | V_T (int p, Vector s, Vector sdsxds, DenseTensor &u) const |
void | V_T (int p, Vector s, Vector sdsxds, real_t dsdsxds, DenseTensor &u, DenseMatrix &du) const |
void | VT_T (int p, Vector s, Vector sds, Vector sdsxds, real_t mu, Vector grad_mu, DenseTensor &u) const |
void | VT_T (int p, Vector s, Vector sds, Vector sdsxds, Vector grad_s2, real_t mu, Vector grad_mu, DenseTensor &u, DenseMatrix &du) const |
void | V_L (int p, Vector sx, const DenseMatrix &grad_sx, Vector sy, const DenseMatrix &grad_sy, real_t t, Vector grad_t, DenseTensor &u) const |
void | V_R (int p, Vector s, const DenseMatrix &grad_s, real_t mu, Vector grad_mu, real_t t, Vector grad_t, DenseMatrix &u) const |
Additional Inherited Members | |
![]() | |
enum | RangeType { UNKNOWN_RANGE_TYPE = -1 , SCALAR , VECTOR } |
Enumeration for range_type and deriv_range_type. More... | |
enum | MapType { UNKNOWN_MAP_TYPE = -1 , VALUE , INTEGRAL , H_DIV , H_CURL } |
Enumeration for MapType: defines how reference functions are mapped to physical space. More... | |
enum | DerivType { NONE , GRAD , DIV , CURL } |
Enumeration for DerivType: defines which derivative method is implemented. More... | |
![]() | |
static bool | IsClosedType (int b_type) |
Return true if the BasisType of b_type is closed (has Quadrature1D points on the boundary). | |
static bool | IsOpenType (int b_type) |
Return true if the BasisType of b_type is open (doesn't have Quadrature1D points on the boundary). | |
static int | VerifyClosed (int b_type) |
Ensure that the BasisType of b_type is closed (has Quadrature1D points on the boundary). | |
static int | VerifyOpen (int b_type) |
Ensure that the BasisType of b_type is open (doesn't have Quadrature1D points on the boundary). | |
static int | VerifyNodal (int b_type) |
Ensure that the BasisType of b_type nodal (satisfies the interpolation property). | |
![]() | |
static bool | CheckZ (real_t z) |
static real_t | lam1 (real_t x, real_t y, real_t z) |
Pyramid "Affine" Coordinates. | |
static real_t | lam2 (real_t x, real_t y, real_t z) |
static real_t | lam3 (real_t x, real_t y, real_t z) |
static real_t | lam4 (real_t x, real_t y, real_t z) |
static real_t | lam5 (real_t x, real_t y, real_t z) |
static Vector | grad_lam1 (real_t x, real_t y, real_t z) |
Gradients of the "Affine" Coordinates. | |
static Vector | grad_lam2 (real_t x, real_t y, real_t z) |
static Vector | grad_lam3 (real_t x, real_t y, real_t z) |
static Vector | grad_lam4 (real_t x, real_t y, real_t z) |
static Vector | grad_lam5 (real_t x, real_t y, real_t z) |
static Vector | lam15 (real_t x, real_t y, real_t z) |
Two component vectors associated with edges touching the apex. | |
static Vector | lam25 (real_t x, real_t y, real_t z) |
static Vector | lam35 (real_t x, real_t y, real_t z) |
static Vector | lam45 (real_t x, real_t y, real_t z) |
static DenseMatrix | grad_lam15 (real_t x, real_t y, real_t z) |
Gradients of the above two component vectors. | |
static DenseMatrix | grad_lam25 (real_t x, real_t y, real_t z) |
static DenseMatrix | grad_lam35 (real_t x, real_t y, real_t z) |
static DenseMatrix | grad_lam45 (real_t x, real_t y, real_t z) |
static Vector | lam15_grad_lam15 (real_t x, real_t y, real_t z) |
Computes | |
static Vector | lam25_grad_lam25 (real_t x, real_t y, real_t z) |
static Vector | lam35_grad_lam35 (real_t x, real_t y, real_t z) |
static Vector | lam45_grad_lam45 (real_t x, real_t y, real_t z) |
static Vector | lam125 (real_t x, real_t y, real_t z) |
Three component vectors associated with triangular faces. | |
static Vector | lam235 (real_t x, real_t y, real_t z) |
static Vector | lam345 (real_t x, real_t y, real_t z) |
static Vector | lam435 (real_t x, real_t y, real_t z) |
static Vector | lam415 (real_t x, real_t y, real_t z) |
static Vector | lam145 (real_t x, real_t y, real_t z) |
static Vector | lam125_grad_lam125 (real_t x, real_t y, real_t z) |
static Vector | lam235_grad_lam235 (real_t x, real_t y, real_t z) |
static Vector | lam345_grad_lam345 (real_t x, real_t y, real_t z) |
static Vector | lam435_grad_lam435 (real_t x, real_t y, real_t z) |
static Vector | lam415_grad_lam415 (real_t x, real_t y, real_t z) |
static Vector | lam145_grad_lam145 (real_t x, real_t y, real_t z) |
static real_t | div_lam125_grad_lam125 (real_t x, real_t y, real_t z) |
Divergences of the above "normal" vector functions divided by 3. | |
static real_t | div_lam235_grad_lam235 (real_t x, real_t y, real_t z) |
static real_t | div_lam345_grad_lam345 (real_t x, real_t y, real_t z) |
static real_t | div_lam435_grad_lam435 (real_t x, real_t y, real_t z) |
static real_t | div_lam415_grad_lam415 (real_t x, real_t y, real_t z) |
static real_t | div_lam145_grad_lam145 (real_t x, real_t y, real_t z) |
static real_t | mu0 (real_t z) |
static real_t | mu1 (real_t z) |
static Vector | grad_mu0 (real_t z) |
static Vector | grad_mu1 (real_t z) |
static Vector | mu01 (real_t z) |
static DenseMatrix | grad_mu01 (real_t z) |
static real_t | mu0 (real_t z, const Vector &xy, unsigned int ab) |
static real_t | mu1 (real_t z, const Vector &xy, unsigned int ab) |
static Vector | grad_mu0 (real_t z, const Vector xy, unsigned int ab) |
static Vector | grad_mu1 (real_t z, const Vector xy, unsigned int ab) |
static Vector | mu01 (real_t z, Vector xy, unsigned int ab) |
static DenseMatrix | grad_mu01 (real_t z, Vector xy, unsigned int ab) |
static Vector | mu01_grad_mu01 (real_t z, Vector xy, unsigned int ab) |
static real_t | nu0 (real_t z, Vector xy, unsigned int ab) |
static real_t | nu1 (real_t z, Vector xy, unsigned int ab) |
static real_t | nu2 (real_t z, Vector xy, unsigned int ab) |
static Vector | grad_nu0 (real_t z, const Vector xy, unsigned int ab) |
static Vector | grad_nu1 (real_t z, const Vector xy, unsigned int ab) |
static Vector | grad_nu2 (real_t z, const Vector xy, unsigned int ab) |
static Vector | nu01 (real_t z, Vector xy, unsigned int ab) |
static Vector | nu12 (real_t z, Vector xy, unsigned int ab) |
static Vector | nu012 (real_t z, Vector xy, unsigned int ab) |
static Vector | nu120 (real_t z, Vector xy, unsigned int ab) |
static DenseMatrix | grad_nu01 (real_t z, Vector xy, unsigned int ab) |
static DenseMatrix | grad_nu012 (real_t z, Vector xy, unsigned int ab) |
static DenseMatrix | grad_nu120 (real_t z, Vector xy, unsigned int ab) |
static Vector | nu01_grad_nu01 (real_t z, Vector xy, unsigned int ab) |
static Vector | nu12_grad_nu12 (real_t z, Vector xy, unsigned int ab) |
static Vector | nu012_grad_nu012 (real_t z, Vector xy, unsigned int ab) |
static void | CalcScaledLegendre (int p, real_t x, real_t t, real_t *u) |
Shifted and Scaled Legendre Polynomials. | |
static void | CalcScaledLegendre (int p, real_t x, real_t t, real_t *u, real_t *dudx, real_t *dudt) |
static void | CalcScaledLegendre (int p, real_t x, real_t t, Vector &u) |
static void | CalcScaledLegendre (int p, real_t x, real_t t, Vector &u, Vector &dudx, Vector &dudt) |
static void | CalcIntegratedLegendre (int p, real_t x, real_t t, real_t *u) |
Integrated Legendre Polynomials. | |
static void | CalcIntegratedLegendre (int p, real_t x, real_t t, real_t *u, real_t *dudx, real_t *dudt) |
static void | CalcIntegratedLegendre (int p, real_t x, real_t t, Vector &u) |
static void | CalcIntegratedLegendre (int p, real_t x, real_t t, Vector &u, Vector &dudx, Vector &dudt) |
static void | CalcHomogenizedScaLegendre (int p, real_t s0, real_t s1, real_t *u) |
static void | CalcHomogenizedScaLegendre (int p, real_t s0, real_t s1, real_t *u, real_t *duds0, real_t *duds1) |
static void | CalcHomogenizedScaLegendre (int p, real_t s0, real_t s1, Vector &u) |
static void | CalcHomogenizedScaLegendre (int p, real_t s0, real_t s1, Vector &u, Vector &duds0, Vector &duds1) |
static void | CalcHomogenizedIntLegendre (int p, real_t t0, real_t t1, real_t *u) |
static void | CalcHomogenizedIntLegendre (int p, real_t t0, real_t t1, real_t *u, real_t *dudt0, real_t *dudt1) |
static void | CalcHomogenizedIntLegendre (int p, real_t t0, real_t t1, Vector &u) |
static void | CalcHomogenizedIntLegendre (int p, real_t t0, real_t t1, Vector &u, Vector &dudt0, Vector &dudt1) |
static void | CalcScaledJacobi (int p, real_t alpha, real_t x, real_t t, real_t *u) |
Shifted and Scaled Jacobi Polynomials. | |
static void | CalcScaledJacobi (int p, real_t alpha, real_t x, real_t t, real_t *u, real_t *dudx, real_t *dudt) |
static void | CalcScaledJacobi (int p, real_t alpha, real_t x, real_t t, Vector &u) |
static void | CalcScaledJacobi (int p, real_t alpha, real_t x, real_t t, Vector &u, Vector &dudx, Vector &dudt) |
static void | CalcIntegratedJacobi (int p, real_t alpha, real_t x, real_t t, real_t *u) |
Integrated Jacobi Polynomials. | |
static void | CalcIntegratedJacobi (int p, real_t alpha, real_t x, real_t t, real_t *u, real_t *dudx, real_t *dudt) |
static void | CalcIntegratedJacobi (int p, real_t alpha, real_t x, real_t t, Vector &u) |
static void | CalcIntegratedJacobi (int p, real_t alpha, real_t x, real_t t, Vector &u, Vector &dudx, Vector &dudt) |
static void | CalcHomogenizedScaJacobi (int p, real_t alpha, real_t t0, real_t t1, real_t *u) |
static void | CalcHomogenizedScaJacobi (int p, real_t alpha, real_t t0, real_t t1, real_t *u, real_t *dudt0, real_t *dudt1) |
static void | CalcHomogenizedScaJacobi (int p, real_t alpha, real_t t0, real_t t1, Vector &u) |
static void | CalcHomogenizedScaJacobi (int p, real_t alpha, real_t t0, real_t t1, Vector &u, Vector &dudt0, Vector &dudt1) |
static void | CalcHomogenizedIntJacobi (int p, real_t alpha, real_t t0, real_t t1, real_t *u) |
static void | CalcHomogenizedIntJacobi (int p, real_t alpha, real_t t0, real_t t1, real_t *u, real_t *dudt0, real_t *dudt1) |
static void | CalcHomogenizedIntJacobi (int p, real_t alpha, real_t t0, real_t t1, Vector &u) |
static void | CalcHomogenizedIntJacobi (int p, real_t alpha, real_t t0, real_t t1, Vector &u, Vector &dudt0, Vector &dudt1) |
static void | phi_E (int p, real_t s0, real_t s1, real_t *u) |
static void | phi_E (int p, real_t s0, real_t s1, real_t *u, real_t *duds0, real_t *duds1) |
static void | phi_E (int p, Vector s, Vector &u) |
static void | phi_E (int p, Vector s, Vector &u, DenseMatrix &duds) |
![]() | |
void | ProjectCurl_2D (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const |
![]() | |
static const ScalarFiniteElement & | CheckScalarFE (const FiniteElement &fe) |
![]() | |
Array< int > | lex_ordering |
![]() | |
int | dim |
Dimension of reference space. | |
int | vdim |
Vector dimension of vector-valued basis functions. | |
int | cdim |
Dimension of curl for vector-valued basis functions. | |
Geometry::Type | geom_type |
Geometry::Type of the reference element. | |
int | func_space |
int | range_type |
int | map_type |
int | deriv_type |
int | deriv_range_type |
int | deriv_map_type |
int | dof |
Number of degrees of freedom. | |
int | order |
Order/degree of the shape functions. | |
int | orders [Geometry::MaxDim] |
Anisotropic orders. | |
IntegrationRule | Nodes |
DenseMatrix | vshape |
Array< DofToQuad * > | dof2quad_array |
Container for all DofToQuad objects created by the FiniteElement. | |
![]() | |
static constexpr real_t | one = 1.0 |
static constexpr real_t | zero = 0.0 |
static constexpr real_t | apex_tol = 1e-8 |
Arbitrary order L2 basis functions defined on pyramid-shaped elements
This implementation is closely based on the finite elements described in section 9.4 of the paper "Orientation embedded high order shape functions for the exact sequence elements of all shapes" by Federico Fuentes, Brendan Keith, Leszek Demkowicz, and Sriram Nagaraj, see https://doi.org/10.1016/j.camwa.2015.04.027.
mfem::L2_FuentesPyramidElement::L2_FuentesPyramidElement | ( | const int | p, |
const int | btype = BasisType::GaussLegendre ) |
|
virtual |
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the given point ip.
Each row of the result DenseMatrix dshape contains the derivatives of one shape function. The size (dof x dim) of dshape must be set in advance.
Implements mfem::FiniteElement.
|
virtual |
Evaluate the values of all shape functions of a scalar finite element in reference space at the given point ip.
The size (dof) of the result Vector shape must be set in advance.
Implements mfem::FiniteElement.