|
MFEM
v3.3
Finite element discretization library
|
Abstract class for Finite Elements. More...
#include <fe.hpp>
Public Types | |
| enum | { SCALAR, VECTOR } |
| Enumeration for RangeType and DerivRangeType. More... | |
| enum | { VALUE, INTEGRAL, H_DIV, H_CURL } |
| Enumeration for MapType: defines how reference functions are mapped to physical space. More... | |
| enum | { NONE, GRAD, DIV, CURL } |
| Enumeration for DerivType: defines which derivative method is implemented. More... | |
Public Member Functions | |
| FiniteElement (int D, int G, int Do, int O, int F=FunctionSpace::Pk) | |
| int | GetDim () const |
| Returns the reference space dimension for the finite element. More... | |
| int | 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. More... | |
| int | Space () const |
| Returns the type of space on each element. More... | |
| int | GetRangeType () const |
| int | GetDerivRangeType () const |
| int | GetMapType () const |
| int | GetDerivType () const |
| int | GetDerivMapType () const |
| 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 |
| 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... | |
| 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 |
| virtual void | CalcHessian (const IntegrationPoint &ip, DenseMatrix &h) const |
| virtual void | GetLocalInterpolation (ElementTransformation &Trans, DenseMatrix &I) const |
| virtual void | Project (Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const |
| virtual void | Project (VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const |
| virtual void | ProjectDelta (int vertex, Vector &dofs) const |
| virtual void | Project (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const |
| virtual void | ProjectGrad (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const |
| virtual void | ProjectCurl (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const |
| virtual void | ProjectDiv (const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &div) const |
| virtual | ~FiniteElement () |
Static Public Member Functions | |
| static int | VerifyClosed (int pt_type) |
| static int | VerifyOpen (int pt_type) |
Protected Attributes | |
| int | Dim |
| Dimension of reference space. More... | |
| int | GeomType |
| Geometry::Type of the reference element. More... | |
| int | Dof |
| Number of degrees of freedom. More... | |
| int | Order |
| Order/degree of the shape functions. More... | |
| int | FuncSpace |
| int | RangeType |
| int | MapType |
| int | DerivType |
| int | DerivRangeType |
| int | DerivMapType |
| IntegrationRule | Nodes |
| DenseMatrix | vshape |
| anonymous enum |
| anonymous enum |
Enumeration for MapType: defines how reference functions are mapped to physical space.
A reference function, uh(xh), can be mapped to a function, u(x), on a general physical element in following ways:
VALUE u(x) = uh(xh)
INTEGRAL u(x) = (1/w) * uh(xh)
H_DIV u(x) = (J/w) * uh(xh)
H_CURL u(x) = J^{-t} * uh(xh) (square J)
H_CURL u(x) = J*(J^t*J)^{-1} * uh(xh) (general J)
where
x = T(xh) is the image of the reference point xh ("x hat"),
J = J(xh) is the Jacobian matrix of the transformation T, and
w = w(xh) = / det(J), for square J,
\ det(J^t*J)^{1/2}, for general J,
is the transformation weight factor. | anonymous enum |
Enumeration for DerivType: defines which derivative method is implemented.
Each FiniteElement class implements only 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. |
| mfem::FiniteElement::FiniteElement | ( | int | D, |
| int | 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 |
|
inlinevirtual |
|
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::ND_TriangleElement, mfem::ND_TetrahedronElement, mfem::ND_QuadrilateralElement, mfem::ND_HexahedronElement, mfem::Nedelec1TetFiniteElement, and mfem::Nedelec1HexFiniteElement.
|
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::RT_TetrahedronElement, mfem::RT_TriangleElement, mfem::RT_HexahedronElement, mfem::RT_QuadrilateralElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::RT2QuadFiniteElement, mfem::RT2TriangleFiniteElement, mfem::RT1QuadFiniteElement, mfem::RT1TriangleFiniteElement, mfem::RT0QuadFiniteElement, and mfem::RT0TriangleFiniteElement.
|
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::NURBS3DFiniteElement, mfem::NURBS2DFiniteElement, mfem::NURBS1DFiniteElement, mfem::L2Pos_TetrahedronElement, mfem::L2_TetrahedronElement, mfem::L2Pos_TriangleElement, mfem::L2_TriangleElement, mfem::L2Pos_HexahedronElement, mfem::L2_HexahedronElement, mfem::L2Pos_QuadrilateralElement, mfem::L2_QuadrilateralElement, mfem::L2Pos_SegmentElement, mfem::L2_SegmentElement, mfem::H1Pos_TetrahedronElement, mfem::H1Pos_TriangleElement, mfem::H1_TetrahedronElement, mfem::H1_TriangleElement, mfem::H1Pos_HexahedronElement, mfem::H1Pos_QuadrilateralElement, mfem::H1Pos_SegmentElement, mfem::H1_HexahedronElement, mfem::H1_QuadrilateralElement, mfem::H1_SegmentElement, mfem::RotTriLinearHexFiniteElement, mfem::RefinedTriLinear3DFiniteElement, mfem::RefinedBiLinear2DFiniteElement, mfem::RefinedLinear3DFiniteElement, mfem::RefinedLinear2DFiniteElement, mfem::RefinedLinear1DFiniteElement, mfem::LagrangeHexFiniteElement, 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::P0TriangleFiniteElement, mfem::Cubic3DFiniteElement, mfem::Cubic2DFiniteElement, mfem::Cubic1DFiniteElement, mfem::BiCubic2DFiniteElement, mfem::GaussBiQuad2DFiniteElement, mfem::BiQuadPos2DFiniteElement, mfem::BiQuad2DFiniteElement, mfem::GaussQuad2DFiniteElement, mfem::Quad2DFiniteElement, mfem::QuadPos1DFiniteElement, mfem::Quad1DFiniteElement, mfem::P1OnQuadFiniteElement, mfem::GaussBiLinear2DFiniteElement, mfem::GaussLinear2DFiniteElement, mfem::BiLinear2DFiniteElement, mfem::Linear2DFiniteElement, mfem::Linear1DFiniteElement, and mfem::PointFiniteElement.
|
virtual |
each row of h contains the upper triangular part of the hessian of one shape function; the order in 2D is {u_xx, u_xy, u_yy}
Reimplemented in mfem::Cubic2DFiniteElement, mfem::BiCubic2DFiniteElement, mfem::Quad2DFiniteElement, and mfem::BiLinear2DFiniteElement.
| void mfem::FiniteElement::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.
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.
| void mfem::FiniteElement::CalcPhysDivShape | ( | ElementTransformation & | Trans, |
| Vector & | divshape | ||
| ) | const |
| 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.
| void mfem::FiniteElement::CalcPhysShape | ( | ElementTransformation & | Trans, |
| Vector & | shape | ||
| ) | const |
|
inline |
Equivalent to the CalcVShape() method with the same arguments.
|
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::NURBS3DFiniteElement, mfem::NURBS2DFiniteElement, mfem::NURBS1DFiniteElement, mfem::ND_SegmentElement, mfem::L2Pos_TetrahedronElement, mfem::L2_TetrahedronElement, mfem::L2Pos_TriangleElement, mfem::L2_TriangleElement, mfem::L2Pos_HexahedronElement, mfem::L2_HexahedronElement, mfem::L2Pos_QuadrilateralElement, mfem::L2_QuadrilateralElement, mfem::L2Pos_SegmentElement, mfem::L2_SegmentElement, mfem::H1Pos_TetrahedronElement, mfem::H1Pos_TriangleElement, mfem::H1_TetrahedronElement, mfem::H1_TriangleElement, mfem::H1Pos_HexahedronElement, mfem::H1Pos_QuadrilateralElement, mfem::H1Pos_SegmentElement, mfem::H1_HexahedronElement, mfem::H1_QuadrilateralElement, mfem::H1_SegmentElement, mfem::RotTriLinearHexFiniteElement, mfem::RefinedTriLinear3DFiniteElement, mfem::RefinedBiLinear2DFiniteElement, mfem::RefinedLinear3DFiniteElement, mfem::RefinedLinear2DFiniteElement, mfem::RefinedLinear1DFiniteElement, mfem::LagrangeHexFiniteElement, 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::P0TriangleFiniteElement, mfem::Cubic3DFiniteElement, mfem::Cubic2DFiniteElement, mfem::Cubic1DFiniteElement, mfem::BiCubic2DFiniteElement, mfem::GaussBiQuad2DFiniteElement, mfem::BiQuadPos2DFiniteElement, mfem::BiQuad2DFiniteElement, mfem::GaussQuad2DFiniteElement, mfem::Quad2DFiniteElement, mfem::QuadPos1DFiniteElement, mfem::Quad1DFiniteElement, mfem::P1OnQuadFiniteElement, mfem::GaussBiLinear2DFiniteElement, mfem::GaussLinear2DFiniteElement, mfem::BiLinear2DFiniteElement, mfem::Linear2DFiniteElement, mfem::Linear1DFiniteElement, and mfem::PointFiniteElement.
|
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::ND_SegmentElement, mfem::ND_TriangleElement, mfem::ND_TetrahedronElement, mfem::ND_QuadrilateralElement, mfem::ND_HexahedronElement, mfem::RT_TetrahedronElement, mfem::RT_TriangleElement, mfem::RT_HexahedronElement, mfem::RT_QuadrilateralElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::RT2QuadFiniteElement, mfem::RT2TriangleFiniteElement, mfem::RT1QuadFiniteElement, mfem::RT1TriangleFiniteElement, mfem::RT0QuadFiniteElement, and mfem::RT0TriangleFiniteElement.
|
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::ND_SegmentElement, mfem::ND_TriangleElement, mfem::ND_TetrahedronElement, mfem::ND_QuadrilateralElement, mfem::ND_HexahedronElement, mfem::RT_TetrahedronElement, mfem::RT_TriangleElement, mfem::RT_HexahedronElement, mfem::RT_QuadrilateralElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::RT2QuadFiniteElement, mfem::RT2TriangleFiniteElement, mfem::RT1QuadFiniteElement, mfem::RT1TriangleFiniteElement, mfem::RT0QuadFiniteElement, and mfem::RT0TriangleFiniteElement.
|
inline |
|
inline |
|
virtual |
Reimplemented in mfem::Linear3DFiniteElement.
|
inline |
Returns the Geometry::Type of the reference element.
|
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::ND_SegmentElement, mfem::ND_TriangleElement, mfem::ND_TetrahedronElement, mfem::ND_QuadrilateralElement, mfem::ND_HexahedronElement, mfem::RT_TetrahedronElement, mfem::RT_TriangleElement, mfem::RT_HexahedronElement, mfem::RT_QuadrilateralElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::RT2QuadFiniteElement, mfem::RT1QuadFiniteElement, mfem::RT1TriangleFiniteElement, mfem::RT0QuadFiniteElement, mfem::RT0TriangleFiniteElement, mfem::BiQuadPos2DFiniteElement, mfem::PositiveFiniteElement, and mfem::NodalFiniteElement.
|
inline |
|
inline |
|
virtual |
Given a coefficient and a transformation, compute its projection (approximation) in the local finite dimensional space in terms of the degrees of freedom.
Reimplemented in mfem::BiQuadPos2DFiniteElement, mfem::PositiveFiniteElement, and mfem::NodalFiniteElement.
|
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)
Reimplemented in mfem::ND_SegmentElement, mfem::ND_TriangleElement, mfem::ND_TetrahedronElement, mfem::ND_QuadrilateralElement, mfem::ND_HexahedronElement, mfem::RT_TetrahedronElement, mfem::RT_TriangleElement, mfem::RT_HexahedronElement, mfem::RT_QuadrilateralElement, mfem::RT0TetFiniteElement, mfem::RT1HexFiniteElement, mfem::RT0HexFiniteElement, mfem::Nedelec1TetFiniteElement, mfem::Nedelec1HexFiniteElement, mfem::RT2QuadFiniteElement, mfem::RT1QuadFiniteElement, mfem::RT1TriangleFiniteElement, mfem::RT0QuadFiniteElement, mfem::RT0TriangleFiniteElement, mfem::BiQuadPos2DFiniteElement, and mfem::NodalFiniteElement.
|
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::ND_SegmentElement, mfem::ND_TriangleElement, mfem::ND_TetrahedronElement, mfem::ND_QuadrilateralElement, mfem::ND_HexahedronElement, mfem::RT_TetrahedronElement, mfem::RT_TriangleElement, mfem::RT_HexahedronElement, mfem::RT_QuadrilateralElement, mfem::PositiveFiniteElement, and mfem::NodalFiniteElement.
|
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::ND_TetrahedronElement, mfem::ND_HexahedronElement, mfem::RT_TetrahedronElement, mfem::RT_TriangleElement, mfem::RT_HexahedronElement, mfem::RT_QuadrilateralElement, mfem::L2_TriangleElement, and mfem::L2_QuadrilateralElement.
|
virtual |
Compute a representation (up to multiplicative constant) for the delta function at the vertex with the given index.
Reimplemented in mfem::L2Pos_TetrahedronElement, mfem::L2_TetrahedronElement, mfem::L2Pos_TriangleElement, mfem::L2_TriangleElement, mfem::L2Pos_HexahedronElement, mfem::L2_HexahedronElement, mfem::L2Pos_QuadrilateralElement, mfem::L2_QuadrilateralElement, mfem::L2Pos_SegmentElement, mfem::L2_SegmentElement, mfem::H1Pos_HexahedronElement, mfem::H1Pos_QuadrilateralElement, mfem::H1Pos_SegmentElement, mfem::H1_HexahedronElement, mfem::H1_QuadrilateralElement, mfem::H1_SegmentElement, mfem::P0HexFiniteElement, mfem::P0TetFiniteElement, mfem::CrouzeixRaviartFiniteElement, mfem::TriLinear3DFiniteElement, mfem::Linear3DFiniteElement, mfem::P0QuadFiniteElement, mfem::P0TriangleFiniteElement, mfem::BiQuadPos2DFiniteElement, mfem::BiQuad2DFiniteElement, mfem::Quad2DFiniteElement, mfem::P1OnQuadFiniteElement, mfem::GaussBiLinear2DFiniteElement, mfem::GaussLinear2DFiniteElement, mfem::BiLinear2DFiniteElement, and mfem::Linear2DFiniteElement.
|
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.
|
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::ND_SegmentElement, mfem::ND_TriangleElement, mfem::ND_TetrahedronElement, mfem::ND_QuadrilateralElement, mfem::ND_HexahedronElement, mfem::RT_TriangleElement, mfem::RT_QuadrilateralElement, and mfem::NodalFiniteElement.
|
inline |
|
inlinestatic |
|
inlinestatic |
|
protected |
|
protected |
|
protected |
Geometry::Type of the reference element.
|
protected |
|
protected |
|
mutableprotected |
1.8.5