MFEM
v4.6.0
Finite element discretization library
|
Abstract class for all finite elements. More...
#include <fe_base.hpp>
Public Types | |
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... | |
Public Member Functions | |
FiniteElement (int D, Geometry::Type G, int Do, int O, int F=FunctionSpace::Pk) | |
Construct FiniteElement with given. More... | |
int | GetDim () const |
Returns the reference space dimension for the finite element. More... | |
int | GetVDim () const |
Returns the vector dimension for vector-valued finite elements. More... | |
int | GetCurlDim () const |
Returns the dimension of the curl for vector-valued finite elements. More... | |
Geometry::Type | GetGeomType () const |
Returns the Geometry::Type of the reference element. More... | |
int | GetDof () const |
Returns the number of degrees of freedom in the finite element. More... | |
int | GetOrder () const |
Returns the order of the finite element. In the case of anisotropic orders, returns the maximum order. More... | |
bool | HasAnisotropicOrders () const |
Returns true if the FiniteElement basis may be using different orders/degrees in different spatial directions. More... | |
const int * | GetAnisotropicOrders () const |
Returns an array containing the anisotropic orders/degrees. More... | |
int | Space () const |
Returns the type of FunctionSpace on the element. More... | |
int | GetRangeType () const |
Returns the FiniteElement::RangeType of the element, one of {SCALAR, VECTOR}. More... | |
int | GetDerivRangeType () const |
Returns the FiniteElement::RangeType of the element derivative, either SCALAR or VECTOR. More... | |
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}. More... | |
int | GetDerivType () const |
Returns the FiniteElement::DerivType of the element describing the spatial derivative method implemented, one of {NONE, GRAD, DIV, CURL}. More... | |
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}. More... | |
virtual void | CalcShape (const IntegrationPoint &ip, Vector &shape) const =0 |
Evaluate the values of all shape functions of a scalar finite element in reference space at the given point ip. More... | |
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. More... | |
virtual void | CalcDShape (const IntegrationPoint &ip, DenseMatrix &dshape) const =0 |
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the given point ip. More... | |
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. More... | |
const IntegrationRule & | GetNodes () const |
Get a const reference to the nodes of the element. More... | |
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. More... | |
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. More... | |
void | CalcPhysVShape (ElementTransformation &Trans, DenseMatrix &shape) const |
Equivalent to the CalcVShape() method with the same arguments. More... | |
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. More... | |
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. More... | |
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. More... | |
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. More... | |
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. More... | |
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. More... | |
virtual void | CalcPhysHessian (ElementTransformation &Trans, DenseMatrix &Hessian) const |
Evaluate the Hessian of all shape functions of a scalar finite element in reference space at the given point ip. More... | |
virtual void | CalcPhysLaplacian (ElementTransformation &Trans, Vector &Laplacian) const |
Evaluate the Laplacian of all shape functions of a scalar finite element in reference space at the given point ip. More... | |
virtual void | CalcPhysLinLaplacian (ElementTransformation &Trans, Vector &Laplacian) const |
virtual void | GetLocalInterpolation (ElementTransformation &Trans, DenseMatrix &I) const |
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base geometry under the given transformation. More... | |
virtual void | GetLocalRestriction (ElementTransformation &Trans, DenseMatrix &R) const |
Return a local restriction matrix R (Dof x Dof) mapping fine dofs to coarse dofs. More... | |
virtual void | GetTransferMatrix (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const |
Return interpolation matrix, I, which maps dofs from a coarse element, fe, to the fine dofs on this finite element. More... | |
virtual void | Project (Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const |
Given a coefficient and a transformation, compute its projection (approximation) in the local finite dimensional space in terms of the degrees of freedom. More... | |
virtual void | Project (VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const |
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) More... | |
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. More... | |
virtual void | ProjectMatrixCoefficient (MatrixCoefficient &mc, ElementTransformation &T, Vector &dofs) const |
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. More... | |
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. More... | |
virtual void | Project (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const |
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. More... | |
virtual void | ProjectGrad (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const |
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. More... | |
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. More... | |
virtual void | ProjectDiv (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &div) const |
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. More... | |
virtual const DofToQuad & | GetDofToQuad (const IntegrationRule &ir, DofToQuad::Mode mode) const |
Return a DofToQuad structure corresponding to the given IntegrationRule using the given DofToQuad::Mode. More... | |
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. More... | |
virtual StatelessDofTransformation * | GetDofTransformation () const |
Return a DoF transformation object for this particular type of basis. More... | |
virtual | ~FiniteElement () |
Deconstruct the FiniteElement. More... | |
Static Public Member Functions | |
static bool | IsClosedType (int b_type) |
Return true if the BasisType of b_type is closed (has Quadrature1D points on the boundary). More... | |
static bool | IsOpenType (int b_type) |
Return true if the BasisType of b_type is open (doesn't have Quadrature1D points on the boundary). More... | |
static int | VerifyClosed (int b_type) |
Ensure that the BasisType of b_type is closed (has Quadrature1D points on the boundary). More... | |
static int | VerifyOpen (int b_type) |
Ensure that the BasisType of b_type is open (doesn't have Quadrature1D points on the boundary). More... | |
static int | VerifyNodal (int b_type) |
Ensure that the BasisType of b_type nodal (satisfies the interpolation property). More... | |
Protected Attributes | |
int | dim |
Dimension of reference space. More... | |
int | vdim |
Vector dimension of vector-valued basis functions. More... | |
int | cdim |
Dimension of curl for vector-valued basis functions. More... | |
Geometry::Type | geom_type |
Geometry::Type of the reference element. More... | |
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. More... | |
int | order |
Order/degree of the shape functions. More... | |
int | orders [Geometry::MaxDim] |
Anisotropic orders. More... | |
IntegrationRule | Nodes |
DenseMatrix | vshape |
Array< DofToQuad * > | dof2quad_array |
Container for all DofToQuad objects created by the FiniteElement. More... | |
Abstract class for all finite elements.
Definition at line 233 of file fe_base.hpp.
Enumeration for DerivType: defines which derivative method is implemented.
Each FiniteElement class implements up to one type of derivative. The value returned by GetDerivType() indicates which derivative method is implemented.
Enumerator | |
---|---|
NONE | No derivatives implemented. |
GRAD | Implements CalcDShape methods. |
DIV | Implements CalcDivShape methods. |
CURL | Implements CalcCurlShape methods. |
Definition at line 292 of file fe_base.hpp.
Enumeration for MapType: defines how reference functions are mapped to physical space.
A reference function \( \hat u(\hat x) \) can be mapped to a function \( u(x) \) on a general physical element in following ways:
Definition at line 269 of file fe_base.hpp.
Enumeration for range_type and deriv_range_type.
Enumerator | |
---|---|
UNKNOWN_RANGE_TYPE | |
SCALAR | |
VECTOR |
Definition at line 257 of file fe_base.hpp.
mfem::FiniteElement::FiniteElement | ( | int | D, |
Geometry::Type | G, | ||
int | Do, | ||
int | O, | ||
int | F = FunctionSpace::Pk |
||
) |
Construct FiniteElement with given.
D | Reference space dimension |
G | Geometry type (of type Geometry::Type) |
Do | Number of degrees of freedom in the FiniteElement |
O | Order/degree of the FiniteElement |
F | FunctionSpace type of the FiniteElement |
Definition at line 23 of file fe_base.cpp.
|
virtual |
Deconstruct the FiniteElement.
Definition at line 500 of file fe_base.cpp.
|
virtual |
Evaluate the curl of all shape functions of a vector finite element in reference space at the given point ip.
Each row of the result DenseMatrix curl_shape contains the components of the curl of one vector shape function. The size (dof x CDim) of curl_shape must be set in advance, where CDim = 3 for dim = 3 and CDim = 1 for dim = 2.
Reimplemented in mfem::Nedelec1PyrFiniteElement, mfem::Nedelec1WdgFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::ND_R2D_QuadrilateralElement, mfem::ND_R2D_TriangleElement, mfem::ND_R2D_SegmentElement, mfem::ND_R1D_SegmentElement, mfem::ND_WedgeElement, mfem::ND_TriangleElement, mfem::ND_TetrahedronElement, mfem::ND_QuadrilateralElement, and mfem::ND_HexahedronElement.
Definition at line 65 of file fe_base.cpp.
|
virtual |
Evaluate the divergence of all shape functions of a vector finite element in reference space at the given point ip.
The size (dof) of the result Vector divshape must be set in advance.
Reimplemented in mfem::RT0PyrFiniteElement, mfem::RT0WdgFiniteElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::RT2QuadFiniteElement, mfem::RT2TriangleFiniteElement, mfem::RT1QuadFiniteElement, mfem::RT1TriangleFiniteElement, mfem::RT0QuadFiniteElement, mfem::RT_R2D_QuadrilateralElement, mfem::RT0TriangleFiniteElement, mfem::RT_R2D_TriangleElement, mfem::RT_R2D_SegmentElement, mfem::RT_R1D_SegmentElement, mfem::RT_WedgeElement, mfem::RT_TetrahedronElement, mfem::RT_TriangleElement, mfem::RT_HexahedronElement, and mfem::RT_QuadrilateralElement.
Definition at line 52 of file fe_base.cpp.
|
pure 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.
Implemented in mfem::RotTriLinearHexFiniteElement, mfem::RefinedTriLinear3DFiniteElement, mfem::RefinedBiLinear2DFiniteElement, mfem::RefinedLinear3DFiniteElement, mfem::RefinedLinear2DFiniteElement, mfem::RefinedLinear1DFiniteElement, mfem::LagrangeHexFiniteElement, mfem::P0PyrFiniteElement, mfem::P0WdgFiniteElement, mfem::P0HexFiniteElement, mfem::P0TetFiniteElement, mfem::P1TetNonConfFiniteElement, mfem::Lagrange1DFiniteElement, mfem::P2SegmentFiniteElement, mfem::P1SegmentFiniteElement, mfem::P0SegmentFiniteElement, mfem::CrouzeixRaviartQuadFiniteElement, mfem::CrouzeixRaviartFiniteElement, mfem::TriLinear3DFiniteElement, mfem::Quadratic3DFiniteElement, mfem::Linear3DFiniteElement, mfem::P0QuadFiniteElement, mfem::L2Pos_WedgeElement, mfem::P0TriangleFiniteElement, mfem::L2Pos_TetrahedronElement, mfem::LinearPyramidFiniteElement, mfem::L2Pos_TriangleElement, mfem::LinearWedgeFiniteElement, mfem::L2Pos_HexahedronElement, mfem::Cubic3DFiniteElement, mfem::L2Pos_QuadrilateralElement, mfem::Cubic2DFiniteElement, mfem::L2Pos_SegmentElement, mfem::Cubic1DFiniteElement, mfem::H1Pos_WedgeElement, mfem::BiCubic2DFiniteElement, mfem::GaussBiQuad2DFiniteElement, mfem::H1Pos_TetrahedronElement, mfem::BiQuad2DFiniteElement, mfem::H1Pos_TriangleElement, mfem::GaussQuad2DFiniteElement, mfem::L2_WedgeElement, mfem::Quad2DFiniteElement, mfem::H1Pos_HexahedronElement, mfem::Quad1DFiniteElement, mfem::H1Pos_QuadrilateralElement, mfem::L2_TetrahedronElement, mfem::NURBS3DFiniteElement, mfem::H1_WedgeElement, mfem::P1OnQuadFiniteElement, mfem::H1Pos_SegmentElement, mfem::GaussBiLinear2DFiniteElement, mfem::L2_TriangleElement, mfem::H1_TetrahedronElement, mfem::QuadPos1DFiniteElement, mfem::NURBS2DFiniteElement, mfem::GaussLinear2DFiniteElement, mfem::H1_TriangleElement, mfem::BiLinear2DFiniteElement, mfem::BiQuadPos2DFiniteElement, mfem::L2_HexahedronElement, mfem::NURBS1DFiniteElement, mfem::H1_HexahedronElement, mfem::Linear2DFiniteElement, mfem::L2_QuadrilateralElement, mfem::H1_QuadrilateralElement, mfem::Linear1DFiniteElement, mfem::H1_SegmentElement, mfem::L2_SegmentElement, mfem::PointFiniteElement, and mfem::H1Ser_QuadrilateralElement.
|
virtual |
Evaluate the Hessians of all shape functions of a scalar finite element in reference space at the given point ip.
Each row of the result DenseMatrix Hessian contains upper triangular part of the Hessian of one shape function. The order in 2D is {u_xx, u_xy, u_yy}. The size (dof x (dim (dim+1)/2) of Hessian must be set in advance.
Reimplemented in mfem::Cubic2DFiniteElement, mfem::BiCubic2DFiniteElement, mfem::Quad2DFiniteElement, mfem::NURBS3DFiniteElement, mfem::H1_TetrahedronElement, mfem::NURBS2DFiniteElement, mfem::H1_TriangleElement, mfem::BiLinear2DFiniteElement, mfem::NURBS1DFiniteElement, mfem::H1_HexahedronElement, mfem::H1_QuadrilateralElement, and mfem::H1_SegmentElement.
Definition at line 101 of file fe_base.cpp.
|
virtual |
Evaluate the curl of all shape functions of a vector finite element in physical space at the point described by Trans.
Each row of the result DenseMatrix curl_shape contains the components of the curl of one vector shape function. The size (dof x CDim) of curl_shape must be set in advance, where CDim = 3 for dim = 3 and CDim = 1 for dim = 2.
Reimplemented in mfem::ND_R2D_FiniteElement, and mfem::ND_R1D_SegmentElement.
Definition at line 71 of file fe_base.cpp.
void mfem::FiniteElement::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.
The size (dof) of the result Vector divshape must be set in advance.
Definition at line 58 of file fe_base.cpp.
void mfem::FiniteElement::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.
Each row of the result DenseMatrix dshape contains the derivatives of one shape function. The size (dof x SDim) of dshape must be set in advance, where SDim >= dim is the physical space dimension as described by Trans.
Definition at line 192 of file fe_base.cpp.
|
virtual |
Evaluate the Hessian of all shape functions of a scalar finite element in reference space at the given point ip.
The size (dof, dim*(dim+1)/2) of Hessian must be set in advance.
Definition at line 288 of file fe_base.cpp.
|
virtual |
Evaluate the Laplacian of all shape functions of a scalar finite element in reference space at the given point ip.
The size (dof) of Laplacian must be set in advance.
Definition at line 203 of file fe_base.cpp.
|
virtual |
Definition at line 244 of file fe_base.cpp.
void mfem::FiniteElement::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.
The size (dof) of the result Vector shape must be set in advance.
Definition at line 182 of file fe_base.cpp.
|
inline |
Equivalent to the CalcVShape() method with the same arguments.
Definition at line 411 of file fe_base.hpp.
|
pure 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.
Implemented in mfem::RotTriLinearHexFiniteElement, mfem::RefinedTriLinear3DFiniteElement, mfem::RefinedBiLinear2DFiniteElement, mfem::RefinedLinear3DFiniteElement, mfem::RefinedLinear2DFiniteElement, mfem::RefinedLinear1DFiniteElement, mfem::LagrangeHexFiniteElement, mfem::P0PyrFiniteElement, mfem::P0WdgFiniteElement, mfem::P0HexFiniteElement, mfem::P0TetFiniteElement, mfem::P1TetNonConfFiniteElement, mfem::Lagrange1DFiniteElement, mfem::P2SegmentFiniteElement, mfem::P1SegmentFiniteElement, mfem::P0SegmentFiniteElement, mfem::CrouzeixRaviartQuadFiniteElement, mfem::CrouzeixRaviartFiniteElement, mfem::TriLinear3DFiniteElement, mfem::Quadratic3DFiniteElement, mfem::Linear3DFiniteElement, mfem::P0QuadFiniteElement, mfem::L2Pos_WedgeElement, mfem::P0TriangleFiniteElement, mfem::L2Pos_TetrahedronElement, mfem::LinearPyramidFiniteElement, mfem::L2Pos_TriangleElement, mfem::LinearWedgeFiniteElement, mfem::L2Pos_HexahedronElement, mfem::ND_SegmentElement, mfem::Cubic3DFiniteElement, mfem::L2Pos_QuadrilateralElement, mfem::Cubic2DFiniteElement, mfem::L2Pos_SegmentElement, mfem::Cubic1DFiniteElement, mfem::H1Pos_WedgeElement, mfem::BiCubic2DFiniteElement, mfem::GaussBiQuad2DFiniteElement, mfem::H1Pos_TetrahedronElement, mfem::BiQuad2DFiniteElement, mfem::H1Pos_TriangleElement, mfem::GaussQuad2DFiniteElement, mfem::L2_WedgeElement, mfem::Quad2DFiniteElement, mfem::H1Pos_HexahedronElement, mfem::H1Pos_QuadrilateralElement, mfem::Quad1DFiniteElement, mfem::L2_TetrahedronElement, mfem::NURBS3DFiniteElement, mfem::H1_WedgeElement, mfem::P1OnQuadFiniteElement, mfem::H1Pos_SegmentElement, mfem::GaussBiLinear2DFiniteElement, mfem::L2_TriangleElement, mfem::H1_TetrahedronElement, mfem::QuadPos1DFiniteElement, mfem::NURBS2DFiniteElement, mfem::GaussLinear2DFiniteElement, mfem::H1_TriangleElement, mfem::BiQuadPos2DFiniteElement, mfem::L2_HexahedronElement, mfem::BiLinear2DFiniteElement, mfem::NURBS1DFiniteElement, mfem::H1_HexahedronElement, mfem::Linear2DFiniteElement, mfem::L2_QuadrilateralElement, mfem::H1_QuadrilateralElement, mfem::Linear1DFiniteElement, mfem::H1_SegmentElement, mfem::L2_SegmentElement, mfem::PointFiniteElement, and mfem::H1Ser_QuadrilateralElement.
|
virtual |
Evaluate the values of all shape functions of a vector finite element in reference space at the given point ip.
Each row of the result DenseMatrix shape contains the components of one vector shape function. The size (dof x dim) of shape must be set in advance.
Reimplemented in mfem::RT0PyrFiniteElement, mfem::RT0WdgFiniteElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::Nedelec1PyrFiniteElement, mfem::Nedelec1WdgFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::ND_R2D_QuadrilateralElement, mfem::RT2QuadFiniteElement, mfem::ND_R2D_TriangleElement, mfem::RT2TriangleFiniteElement, mfem::RT1QuadFiniteElement, mfem::RT1TriangleFiniteElement, mfem::ND_R2D_SegmentElement, mfem::RT0QuadFiniteElement, mfem::RT_R2D_QuadrilateralElement, mfem::RT0TriangleFiniteElement, mfem::RT_R2D_TriangleElement, mfem::ND_R1D_SegmentElement, mfem::ND_R1D_PointElement, mfem::RT_R2D_SegmentElement, mfem::ND_WedgeElement, mfem::RT_R1D_SegmentElement, mfem::ND_SegmentElement, mfem::RT_WedgeElement, mfem::ND_TriangleElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::ND_QuadrilateralElement, mfem::RT_HexahedronElement, mfem::RT_QuadrilateralElement, and mfem::ND_HexahedronElement.
Definition at line 40 of file fe_base.cpp.
|
virtual |
Evaluate the values of all shape functions of a vector finite element in physical space at the point described by Trans.
Each row of the result DenseMatrix shape contains the components of one vector shape function. The size (dof x SDim) of shape must be set in advance, where SDim >= dim is the physical space dimension as described by Trans.
Reimplemented in mfem::RT0PyrFiniteElement, mfem::RT0WdgFiniteElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::Nedelec1PyrFiniteElement, mfem::Nedelec1WdgFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::RT2QuadFiniteElement, mfem::RT2TriangleFiniteElement, mfem::RT1QuadFiniteElement, mfem::ND_R2D_FiniteElement, mfem::RT1TriangleFiniteElement, mfem::ND_R2D_SegmentElement, mfem::RT0QuadFiniteElement, mfem::RT0TriangleFiniteElement, mfem::ND_R1D_SegmentElement, mfem::RT_R2D_FiniteElement, mfem::ND_R1D_PointElement, mfem::RT_R2D_SegmentElement, mfem::ND_WedgeElement, mfem::RT_R1D_SegmentElement, mfem::ND_SegmentElement, mfem::RT_WedgeElement, mfem::ND_TriangleElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::ND_QuadrilateralElement, mfem::RT_HexahedronElement, mfem::RT_QuadrilateralElement, and mfem::ND_HexahedronElement.
Definition at line 46 of file fe_base.cpp.
|
inline |
Returns an array containing the anisotropic orders/degrees.
Definition at line 334 of file fe_base.hpp.
|
inline |
Returns the dimension of the curl for vector-valued finite elements.
Definition at line 317 of file fe_base.hpp.
|
inline |
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}.
Definition at line 359 of file fe_base.hpp.
|
inline |
Returns the FiniteElement::RangeType of the element derivative, either SCALAR or VECTOR.
Definition at line 344 of file fe_base.hpp.
|
inline |
Returns the FiniteElement::DerivType of the element describing the spatial derivative method implemented, one of {NONE, GRAD, DIV, CURL}.
Definition at line 354 of file fe_base.hpp.
|
inline |
Returns the reference space dimension for the finite element.
Definition at line 311 of file fe_base.hpp.
|
inline |
Returns the number of degrees of freedom in the finite element.
Definition at line 323 of file fe_base.hpp.
|
virtual |
Return a DofToQuad structure corresponding to the given IntegrationRule using the given DofToQuad::Mode.
See the documentation for DofToQuad for more details.
Reimplemented in mfem::VectorTensorFiniteElement, mfem::NodalTensorFiniteElement, and mfem::PositiveTensorFiniteElement.
Definition at line 365 of file fe_base.cpp.
|
inlinevirtual |
Return a DoF transformation object for this particular type of basis.
Reimplemented in mfem::ND_WedgeElement, mfem::ND_TriangleElement, and mfem::ND_TetrahedronElement.
Definition at line 598 of file fe_base.hpp.
|
virtual |
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.
Reimplemented in mfem::Linear3DFiniteElement, mfem::LinearPyramidFiniteElement, and mfem::LinearWedgeFiniteElement.
Definition at line 96 of file fe_base.cpp.
|
virtual |
Return the mapping from lexicographic face DOFs to lexicographic element DOFs for the given local face face_id.
Given the ith DOF (lexicographically ordered) on the face referenced by face_id, face_map[i] gives the corresponding index of the DOF in the element (also lexicographically ordered).
Reimplemented in mfem::NodalTensorFiniteElement, mfem::ND_QuadrilateralElement, mfem::RT_HexahedronElement, mfem::ND_HexahedronElement, mfem::RT_QuadrilateralElement, and mfem::PositiveTensorFiniteElement.
Definition at line 494 of file fe_base.cpp.
|
inline |
Returns the Geometry::Type of the reference element.
Definition at line 320 of file fe_base.hpp.
|
virtual |
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base geometry under the given transformation.
Reimplemented in mfem::RT0PyrFiniteElement, mfem::RT0WdgFiniteElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::Nedelec1PyrFiniteElement, mfem::Nedelec1WdgFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::NodalFiniteElement, mfem::RT2QuadFiniteElement, mfem::RT1QuadFiniteElement, mfem::ND_R2D_FiniteElement, mfem::RT1TriangleFiniteElement, mfem::ND_R2D_SegmentElement, mfem::RT0QuadFiniteElement, mfem::RT0TriangleFiniteElement, mfem::ND_R1D_SegmentElement, mfem::RT_R2D_FiniteElement, mfem::RT_R2D_SegmentElement, mfem::ND_WedgeElement, mfem::ND_SegmentElement, mfem::RT_WedgeElement, mfem::ND_TriangleElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::ND_QuadrilateralElement, mfem::RT_HexahedronElement, mfem::BiQuadPos2DFiniteElement, mfem::ND_HexahedronElement, mfem::RT_QuadrilateralElement, mfem::PositiveFiniteElement, and mfem::H1Ser_QuadrilateralElement.
Definition at line 107 of file fe_base.cpp.
|
virtual |
Return a local restriction matrix R (Dof x Dof) mapping fine dofs to coarse dofs.
The fine element is the image of the base geometry under the given transformation, Trans.
The assumption in this method is that a subset of the coarse dofs can be expressed only in terms of the dofs of the given fine element.
Rows in R corresponding to coarse dofs that cannot be expressed in terms of the fine dofs will be marked as invalid by setting the first entry (column 0) in the row to infinity().
This method assumes that the dimensions of R are set before it is called.
Reimplemented in mfem::NodalFiniteElement, mfem::ND_R2D_FiniteElement, mfem::ND_R2D_SegmentElement, mfem::ND_R1D_SegmentElement, mfem::RT_R2D_FiniteElement, mfem::RT_R2D_SegmentElement, mfem::ND_WedgeElement, mfem::ND_SegmentElement, mfem::RT_WedgeElement, mfem::ND_TriangleElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::L2_TetrahedronElement, mfem::ND_QuadrilateralElement, mfem::L2_TriangleElement, mfem::RT_HexahedronElement, mfem::L2_HexahedronElement, mfem::L2_QuadrilateralElement, mfem::ND_HexahedronElement, mfem::RT_QuadrilateralElement, mfem::PositiveFiniteElement, and mfem::L2_SegmentElement.
Definition at line 113 of file fe_base.cpp.
|
inline |
Returns the FiniteElement::MapType of the element describing how reference functions are mapped to physical space, one of {VALUE, INTEGRAL H_DIV, H_CURL}.
Definition at line 349 of file fe_base.hpp.
|
inline |
Get a const reference to the nodes of the element.
Definition at line 389 of file fe_base.hpp.
|
inline |
Returns the order of the finite element. In the case of anisotropic orders, returns the maximum order.
Definition at line 327 of file fe_base.hpp.
|
inline |
Returns the FiniteElement::RangeType of the element, one of {SCALAR, VECTOR}.
Definition at line 340 of file fe_base.hpp.
|
virtual |
Return interpolation matrix, I, which maps dofs from a coarse element, fe, to the fine dofs on this finite element.
Trans represents the mapping from the reference element of this element into a subset of the reference space of the element fe, thus allowing the "coarse" FiniteElement to be different from the "fine" FiniteElement as when h-refinement is combined with p-refinement or p-derefinement. It is assumed that both finite elements use the same FiniteElement::MapType.
Reimplemented in mfem::NodalTensorFiniteElement, mfem::NodalFiniteElement, mfem::ND_R2D_FiniteElement, mfem::ND_R2D_SegmentElement, mfem::ND_R1D_SegmentElement, mfem::RT_R2D_FiniteElement, mfem::RT_R2D_SegmentElement, mfem::ND_WedgeElement, mfem::ND_SegmentElement, mfem::RT_WedgeElement, mfem::ND_TriangleElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::ND_QuadrilateralElement, mfem::RT_HexahedronElement, mfem::ND_HexahedronElement, mfem::RT_QuadrilateralElement, and mfem::PositiveFiniteElement.
Definition at line 119 of file fe_base.cpp.
|
inline |
Returns the vector dimension for vector-valued finite elements.
Definition at line 314 of file fe_base.hpp.
|
inline |
Returns true if the FiniteElement basis may be using different orders/degrees in different spatial directions.
Definition at line 331 of file fe_base.hpp.
|
inlinestatic |
Return true if the BasisType of b_type is closed (has Quadrature1D points on the boundary).
Definition at line 606 of file fe_base.hpp.
|
inlinestatic |
Return true if the BasisType of b_type is open (doesn't have Quadrature1D points on the boundary).
Definition at line 615 of file fe_base.hpp.
|
virtual |
Given a coefficient and a transformation, compute its projection (approximation) in the local finite dimensional space in terms of the degrees of freedom.
The approximation used to project is usually local interpolation of degrees of freedom. The derived class could use other methods not implemented yet, e.g. local L2 projection.
Reimplemented in mfem::NodalFiniteElement, mfem::L2_HexahedronElement, mfem::BiQuadPos2DFiniteElement, mfem::L2_QuadrilateralElement, and mfem::PositiveFiniteElement.
Definition at line 126 of file fe_base.cpp.
|
virtual |
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)
The approximation used to project is usually local interpolation of degrees of freedom. The derived class could use other methods not implemented yet, e.g. local L2 projection.
Reimplemented in mfem::RT0PyrFiniteElement, mfem::RT0WdgFiniteElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::Nedelec1PyrFiniteElement, mfem::Nedelec1WdgFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::NodalFiniteElement, mfem::RT2QuadFiniteElement, mfem::ND_R2D_FiniteElement, mfem::RT1QuadFiniteElement, mfem::RT1TriangleFiniteElement, mfem::ND_R2D_SegmentElement, mfem::RT0QuadFiniteElement, mfem::RT0TriangleFiniteElement, mfem::ND_R1D_SegmentElement, mfem::RT_R2D_FiniteElement, mfem::ND_WedgeElement, mfem::RT_R1D_SegmentElement, mfem::ND_SegmentElement, mfem::RT_WedgeElement, mfem::ND_TriangleElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::ND_QuadrilateralElement, mfem::RT_HexahedronElement, mfem::BiQuadPos2DFiniteElement, mfem::ND_HexahedronElement, mfem::RT_QuadrilateralElement, and mfem::PositiveFiniteElement.
Definition at line 132 of file fe_base.cpp.
|
virtual |
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.
Reimplemented in mfem::NodalFiniteElement, mfem::ND_R2D_FiniteElement, mfem::ND_R1D_SegmentElement, mfem::RT_R2D_FiniteElement, mfem::ND_WedgeElement, mfem::RT_R1D_SegmentElement, mfem::ND_SegmentElement, mfem::RT_WedgeElement, mfem::ND_TriangleElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::ND_QuadrilateralElement, mfem::RT_HexahedronElement, mfem::ND_HexahedronElement, mfem::RT_QuadrilateralElement, and mfem::PositiveFiniteElement.
Definition at line 155 of file fe_base.cpp.
|
virtual |
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.
Reimplemented in mfem::RT0PyrFiniteElement, mfem::RT0WdgFiniteElement, mfem::ND_R1D_SegmentElement, mfem::RT_R2D_FiniteElement, mfem::ND_WedgeElement, mfem::RT_R1D_SegmentElement, mfem::RT_WedgeElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::RT_HexahedronElement, mfem::L2_TriangleElement, mfem::ND_HexahedronElement, mfem::RT_QuadrilateralElement, and mfem::L2_QuadrilateralElement.
Definition at line 168 of file fe_base.cpp.
|
virtual |
Project a delta function centered on the given vertex in the local finite dimensional space represented by the dofs.
Reimplemented in mfem::P0PyrFiniteElement, mfem::P0WdgFiniteElement, mfem::P0HexFiniteElement, mfem::P0TetFiniteElement, mfem::CrouzeixRaviartFiniteElement, mfem::TriLinear3DFiniteElement, mfem::Linear3DFiniteElement, mfem::P0QuadFiniteElement, mfem::P0TriangleFiniteElement, mfem::L2Pos_TetrahedronElement, mfem::LinearPyramidFiniteElement, mfem::L2Pos_TriangleElement, mfem::LinearWedgeFiniteElement, mfem::L2Pos_HexahedronElement, mfem::L2Pos_QuadrilateralElement, mfem::L2Pos_SegmentElement, mfem::BiQuad2DFiniteElement, mfem::Quad2DFiniteElement, mfem::H1Pos_HexahedronElement, mfem::H1Pos_QuadrilateralElement, mfem::L2_TetrahedronElement, mfem::P1OnQuadFiniteElement, mfem::H1Pos_SegmentElement, mfem::GaussBiLinear2DFiniteElement, mfem::L2_TriangleElement, mfem::GaussLinear2DFiniteElement, mfem::BiQuadPos2DFiniteElement, mfem::BiLinear2DFiniteElement, mfem::L2_HexahedronElement, mfem::H1_HexahedronElement, mfem::Linear2DFiniteElement, mfem::L2_QuadrilateralElement, mfem::H1_QuadrilateralElement, mfem::H1_SegmentElement, and mfem::L2_SegmentElement.
Definition at line 150 of file fe_base.cpp.
|
virtual |
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.
Reimplemented in mfem::NodalFiniteElement, mfem::L2_HexahedronElement, and mfem::L2_QuadrilateralElement.
Definition at line 175 of file fe_base.cpp.
|
virtual |
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.
Reimplemented in mfem::ND_R1D_SegmentElement, mfem::ND_TriangleElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::ND_QuadrilateralElement, mfem::RT_HexahedronElement, mfem::ND_HexahedronElement, and mfem::RT_QuadrilateralElement.
Definition at line 138 of file fe_base.cpp.
|
virtual |
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.
Reimplemented in mfem::Nedelec1PyrFiniteElement, mfem::Nedelec1WdgFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::NodalFiniteElement, mfem::ND_R2D_FiniteElement, mfem::ND_R1D_SegmentElement, mfem::ND_WedgeElement, mfem::ND_SegmentElement, mfem::ND_TriangleElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::ND_QuadrilateralElement, mfem::ND_HexahedronElement, and mfem::RT_QuadrilateralElement.
Definition at line 161 of file fe_base.cpp.
|
virtual |
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.
Reimplemented in mfem::NodalFiniteElement, mfem::ND_R1D_SegmentElement, mfem::ND_WedgeElement, mfem::ND_SegmentElement, mfem::RT_WedgeElement, mfem::ND_TriangleElement, mfem::RT_TetrahedronElement, mfem::ND_TetrahedronElement, mfem::RT_TriangleElement, mfem::ND_QuadrilateralElement, mfem::RT_HexahedronElement, mfem::ND_HexahedronElement, and mfem::RT_QuadrilateralElement.
Definition at line 144 of file fe_base.cpp.
|
inline |
Returns the type of FunctionSpace on the element.
Definition at line 337 of file fe_base.hpp.
|
inlinestatic |
Ensure that the BasisType of b_type is closed (has Quadrature1D points on the boundary).
Definition at line 624 of file fe_base.hpp.
|
inlinestatic |
Ensure that the BasisType of b_type nodal (satisfies the interpolation property).
Definition at line 641 of file fe_base.hpp.
|
inlinestatic |
Ensure that the BasisType of b_type is open (doesn't have Quadrature1D points on the boundary).
Definition at line 633 of file fe_base.hpp.
|
protected |
Dimension of curl for vector-valued basis functions.
Definition at line 238 of file fe_base.hpp.
|
protected |
Definition at line 240 of file fe_base.hpp.
|
protected |
Definition at line 240 of file fe_base.hpp.
|
protected |
Definition at line 240 of file fe_base.hpp.
|
protected |
Dimension of reference space.
Definition at line 236 of file fe_base.hpp.
|
mutableprotected |
Number of degrees of freedom.
Definition at line 243 of file fe_base.hpp.
Container for all DofToQuad objects created by the FiniteElement.
Multiple DofToQuad objects may be needed when different quadrature rules or different DofToQuad::Mode are used.
Definition at line 253 of file fe_base.hpp.
|
protected |
Definition at line 240 of file fe_base.hpp.
|
protected |
Geometry::Type of the reference element.
Definition at line 239 of file fe_base.hpp.
|
protected |
Definition at line 240 of file fe_base.hpp.
|
protected |
Definition at line 246 of file fe_base.hpp.
|
mutableprotected |
Order/degree of the shape functions.
Definition at line 243 of file fe_base.hpp.
|
mutableprotected |
Anisotropic orders.
Definition at line 245 of file fe_base.hpp.
|
protected |
Definition at line 240 of file fe_base.hpp.
|
protected |
Vector dimension of vector-valued basis functions.
Definition at line 237 of file fe_base.hpp.
|
mutableprotected |
Definition at line 248 of file fe_base.hpp.