LehrFEM++ 1.0.0
A simple Finite Element Library for teaching
Loading...
Searching...
No Matches
Public Member Functions | List of all members
lf::uscalfe::FeLagrangeO2Segment< SCALAR > Class Template Reference

Quadratic Lagrangian finite element on a line segment. More...

#include <lf/uscalfe/uscalfe.h>

Inheritance diagram for lf::uscalfe::FeLagrangeO2Segment< SCALAR >:
lf::fe::ScalarReferenceFiniteElement< SCALAR >

Public Member Functions

 FeLagrangeO2Segment (const FeLagrangeO2Segment &)=default
 
 FeLagrangeO2Segment (FeLagrangeO2Segment &&) noexcept=default
 
FeLagrangeO2Segmentoperator= (const FeLagrangeO2Segment &)=default
 
FeLagrangeO2Segmentoperator= (FeLagrangeO2Segment &&) noexcept=default
 
 FeLagrangeO2Segment ()=default
 
 ~FeLagrangeO2Segment () override=default
 
base::RefEl RefEl () const override
 Tells the type of reference cell underlying the parametric finite element.
 
unsigned Degree () const override
 Quadratic Lagrangian finite element rely on polynomials of degree 2.
 
size_type NumRefShapeFunctions () const override
 Three local shape functions are associated with an edge.
 
size_type NumRefShapeFunctions (dim_t codim) const override
 One shape function attached to each node and one interior shape function.
 
size_type NumRefShapeFunctions (dim_t codim, sub_idx_t) const override
 One shape function attached to each node and one interior shape function.
 
Eigen::Matrix< SCALAR, Eigen::Dynamic, Eigen::Dynamic > EvalReferenceShapeFunctions (const Eigen::MatrixXd &refcoords) const override
 Evaluation of shape functions on reference segment (unit interval)
 
Eigen::Matrix< SCALAR, Eigen::Dynamic, Eigen::Dynamic > GradientsReferenceShapeFunctions (const Eigen::MatrixXd &refcoords) const override
 Evaluation of derivatives of local shape function on reference interval.
 
Eigen::MatrixXd EvaluationNodes () const override
 Evaluation nodes are the endpoints of the segment and its midpoint.
 
size_type NumEvaluationNodes () const override
 Three evaluation nodes, same number as local shape functions.
 
- Public Member Functions inherited from lf::fe::ScalarReferenceFiniteElement< SCALAR >
virtual ~ScalarReferenceFiniteElement ()=default
 
dim_t Dimension () const
 Returns the spatial dimension of the reference cell.
 
virtual Eigen::Matrix< SCALAR, 1, Eigen::Dynamic > NodalValuesToDofs (const Eigen::Matrix< SCALAR, 1, Eigen::Dynamic > &nodvals) const
 Computes the linear combination of reference shape functions matching function values at evaluation nodes.
 

Additional Inherited Members

- Public Types inherited from lf::fe::ScalarReferenceFiniteElement< SCALAR >
using Scalar = SCALAR
 The underlying scalar type.
 
- Protected Member Functions inherited from lf::fe::ScalarReferenceFiniteElement< SCALAR >
 ScalarReferenceFiniteElement ()=default
 
 ScalarReferenceFiniteElement (const ScalarReferenceFiniteElement &)=default
 
 ScalarReferenceFiniteElement (ScalarReferenceFiniteElement &&) noexcept=default
 
ScalarReferenceFiniteElementoperator= (const ScalarReferenceFiniteElement &)=default
 
ScalarReferenceFiniteElementoperator= (ScalarReferenceFiniteElement &&) noexcept=default
 

Detailed Description

template<typename SCALAR>
class lf::uscalfe::FeLagrangeO2Segment< SCALAR >

Quadratic Lagrangian finite element on a line segment.

This is a specialization of ScalarReferenceFiniteElement for quadratic Lagrangian finite elements on a segment.

The first two shape functions are associated with the vertices of the segment. The last one is an interior shape function. This complies with the ordering convention of LehreFEM++ for local shape functions, see Paragraph 2.7.4.11.

See also
fe::ScalarReferenceFiniteElement

Definition at line 563 of file lagr_fe.h.

Constructor & Destructor Documentation

◆ FeLagrangeO2Segment() [1/3]

template<typename SCALAR >
lf::uscalfe::FeLagrangeO2Segment< SCALAR >::FeLagrangeO2Segment ( const FeLagrangeO2Segment< SCALAR > & )
default

◆ FeLagrangeO2Segment() [2/3]

template<typename SCALAR >
lf::uscalfe::FeLagrangeO2Segment< SCALAR >::FeLagrangeO2Segment ( FeLagrangeO2Segment< SCALAR > && )
defaultnoexcept

◆ FeLagrangeO2Segment() [3/3]

template<typename SCALAR >
lf::uscalfe::FeLagrangeO2Segment< SCALAR >::FeLagrangeO2Segment ( )
default

◆ ~FeLagrangeO2Segment()

template<typename SCALAR >
lf::uscalfe::FeLagrangeO2Segment< SCALAR >::~FeLagrangeO2Segment ( )
overridedefault

Member Function Documentation

◆ Degree()

template<typename SCALAR >
unsigned lf::uscalfe::FeLagrangeO2Segment< SCALAR >::Degree ( ) const
inlineoverridevirtual

Quadratic Lagrangian finite element rely on polynomials of degree 2.

Implements lf::fe::ScalarReferenceFiniteElement< SCALAR >.

Definition at line 579 of file lagr_fe.h.

◆ EvalReferenceShapeFunctions()

template<typename SCALAR >
Eigen::Matrix< SCALAR, Eigen::Dynamic, Eigen::Dynamic > lf::uscalfe::FeLagrangeO2Segment< SCALAR >::EvalReferenceShapeFunctions ( const Eigen::MatrixXd & refcoords) const
inlineoverridevirtual

Evaluation of shape functions on reference segment (unit interval)

The local shape functions are associated with a segment in the case of quadratic Lagrangian finite elements. In reference coordinates those are:

\[ \hat{b}^1(x) = 2(1-x)(\frac12-x)\;,\quad \hat{b}^2(x) = 2x(x-\frac12)\;,\quad \hat{b}^3(x) = 4(1-x)x\;. \]

This is a cardinal basis of the space of quadratic univariate polynomials with respect to the points \(x=0,\frac12,1\).

See also
fe::ScalarReferenceFiniteElement::EvalReferenceShapeFunctions()

Implements lf::fe::ScalarReferenceFiniteElement< SCALAR >.

Definition at line 627 of file lagr_fe.h.

◆ EvaluationNodes()

template<typename SCALAR >
Eigen::MatrixXd lf::uscalfe::FeLagrangeO2Segment< SCALAR >::EvaluationNodes ( ) const
inlineoverridevirtual

Evaluation nodes are the endpoints of the segment and its midpoint.

Returns reference coordinates of "evaluation nodes" for evaluation of parametric degrees of freedom, nodal interpolation in the simplest case.

Returns
A d x N matrix, where d is the dimension of the reference cell, and N the number of interpolation nodes. The columns of this matrix contain their reference coordinates.

Every parametric scalar finite element implicitly defines a local interpolation operator by duality with the reference shape functions. This interpolation operator can be realized through evaluations at certain evaluation nodes, which are provided by this method.

Unisolvence

The evaluation points must satisfy the following requirement: If the values of a function belonging to the span of the reference shape functions are known in the evaluation nodes, then this function is uniquely determined. This entails that the number of evaluation nodes must be at least as big as the number of reference shape functions.

Note
It is not required that any vector a values at evaluation nodes can be produced by a suitable linear combination of reference shape functions. For instance, this will not be possible, if there are more evaluation points than reference shape functions. If both sets have the same size, however, the interpolation problem has a solution for any vector of values at the evluation points.

Example: Principal lattice

For triangular Lagrangian finite elements of degree p the evaluation nodes, usually called "interpolation nodes" in this context, can be chosen as \(\left(\frac{j}{p},\frac{k}{p}\right),\; 0\leq j,k \leq p, j+k\leq p\).

Moment-based degrees of freedom

For some finite element spaces the interpolation functional may be defined based on integrals over edges. In this case the evaluation nodes will be quadrature nodes for the approximate evaluation of these integrals.

The quadrature rule must be exact for the polynomials contained in the local finite element spaces.

Implements lf::fe::ScalarReferenceFiniteElement< SCALAR >.

Definition at line 673 of file lagr_fe.h.

◆ GradientsReferenceShapeFunctions()

template<typename SCALAR >
Eigen::Matrix< SCALAR, Eigen::Dynamic, Eigen::Dynamic > lf::uscalfe::FeLagrangeO2Segment< SCALAR >::GradientsReferenceShapeFunctions ( const Eigen::MatrixXd & refcoords) const
inlineoverridevirtual

Evaluation of derivatives of local shape function on reference interval.

For the passed points return the derivative of the three local shape functions on the reference interval in the rows of a matrix.

See also
FeLagrangeO2Segment::EvalReferenceShapeFunctions()
fe::ScalarReferenceFiniteElement::GradientsReferenceShapeFunctions()

Implements lf::fe::ScalarReferenceFiniteElement< SCALAR >.

Definition at line 653 of file lagr_fe.h.

◆ NumEvaluationNodes()

template<typename SCALAR >
size_type lf::uscalfe::FeLagrangeO2Segment< SCALAR >::NumEvaluationNodes ( ) const
inlineoverridevirtual

Three evaluation nodes, same number as local shape functions.

Tell the number of evaluation (interpolation) nodes.

Implements lf::fe::ScalarReferenceFiniteElement< SCALAR >.

Definition at line 683 of file lagr_fe.h.

References lf::uscalfe::FeLagrangeO2Segment< SCALAR >::NumRefShapeFunctions().

◆ NumRefShapeFunctions() [1/3]

template<typename SCALAR >
size_type lf::uscalfe::FeLagrangeO2Segment< SCALAR >::NumRefShapeFunctions ( ) const
inlineoverridevirtual

◆ NumRefShapeFunctions() [2/3]

template<typename SCALAR >
size_type lf::uscalfe::FeLagrangeO2Segment< SCALAR >::NumRefShapeFunctions ( dim_t codim) const
inlineoverridevirtual

One shape function attached to each node and one interior shape function.

Hence (sub-)entities of any co-dimension possess exactly one local shape function.

See also
fe::ScalarReferenceFiniteElement::NumRefShapeFunctions(dim_t)

Reimplemented from lf::fe::ScalarReferenceFiniteElement< SCALAR >.

Definition at line 596 of file lagr_fe.h.

◆ NumRefShapeFunctions() [3/3]

template<typename SCALAR >
size_type lf::uscalfe::FeLagrangeO2Segment< SCALAR >::NumRefShapeFunctions ( dim_t codim,
sub_idx_t  ) const
inlineoverridevirtual

One shape function attached to each node and one interior shape function.

See also
fe::ScalarReferenceFiniteElement::NumRefShapeFunctions(dim_t,sub_idx_t)

Implements lf::fe::ScalarReferenceFiniteElement< SCALAR >.

Definition at line 606 of file lagr_fe.h.

◆ operator=() [1/2]

template<typename SCALAR >
FeLagrangeO2Segment & lf::uscalfe::FeLagrangeO2Segment< SCALAR >::operator= ( const FeLagrangeO2Segment< SCALAR > & )
default

◆ operator=() [2/2]

template<typename SCALAR >
FeLagrangeO2Segment & lf::uscalfe::FeLagrangeO2Segment< SCALAR >::operator= ( FeLagrangeO2Segment< SCALAR > && )
defaultnoexcept

◆ RefEl()

template<typename SCALAR >
base::RefEl lf::uscalfe::FeLagrangeO2Segment< SCALAR >::RefEl ( ) const
inlineoverridevirtual

Tells the type of reference cell underlying the parametric finite element.

Implements lf::fe::ScalarReferenceFiniteElement< SCALAR >.

Definition at line 573 of file lagr_fe.h.

References lf::base::RefEl::kSegment().


The documentation for this class was generated from the following file: