Public Member Functions | |
FE_RaviartThomasNodal (const unsigned int p) | |
virtual std::string | get_name () const |
virtual FiniteElement< dim > * | clone () const |
virtual void | interpolate (std::vector< double > &local_dofs, const std::vector< double > &values) const |
virtual void | interpolate (std::vector< double > &local_dofs, const std::vector< Vector< double > > &values, unsigned int offset=0) const |
virtual void | interpolate (std::vector< double > &local_dofs, const VectorSlice< const std::vector< std::vector< double > > > &values) const |
virtual void | get_face_interpolation_matrix (const FiniteElement< dim > &source, FullMatrix< double > &matrix) const |
virtual void | get_subface_interpolation_matrix (const FiniteElement< dim > &source, const unsigned int subface, FullMatrix< double > &matrix) const |
virtual bool | hp_constraints_are_implemented () const |
virtual std::vector< std::pair < unsigned int, unsigned int > > | hp_vertex_dof_identities (const FiniteElement< dim > &fe_other) const |
virtual std::vector< std::pair < unsigned int, unsigned int > > | hp_line_dof_identities (const FiniteElement< dim > &fe_other) const |
virtual std::vector< std::pair < unsigned int, unsigned int > > | hp_quad_dof_identities (const FiniteElement< dim > &fe_other) const |
virtual FiniteElementDomination::Domination | compare_for_face_domination (const FiniteElement< dim > &fe_other) const |
Private Member Functions | |
virtual bool | has_support_on_face (const unsigned int shape_index, const unsigned int face_index) const |
void | initialize_support_points (const unsigned int rt_degree) |
Static Private Member Functions | |
static std::vector< unsigned int > | get_dpo_vector (const unsigned int degree) |
static std::vector< bool > | get_ria_vector (const unsigned int degree) |
For this Raviart-Thomas element, the node values are not cell and face moments with respect to certain polynomials, but the values in quadrature points. Following the general scheme for numbering degrees of freedom, the node values on edges are first, edge by edge, according to the natural ordering of the edges of a cell. The interior degrees of freedom are last.
For an RT-element of degree k, we choose (k+1)d-1 Gauss points on each face. These points are ordered lexicographically with respect to the orientation of the face. This way, the normal component which is in Qk is uniquely determined. Furthermore, since this Gauss-formula is exact on Q2k+1, these node values correspond to the exact integration of the moments of the RT-space.
In the interior of the cells, the moments are with respect to an anisotropic Qk space, where the test functions are one degree lower in the direction corresponding to the vector component under consideration. This is emulated by using an anisotropic Gauss formula for integration.
FE_RaviartThomasNodal< dim >::FE_RaviartThomasNodal | ( | const unsigned int | p | ) |
Constructor for the Raviart-Thomas element of degree p
.
virtual std::string FE_RaviartThomasNodal< dim >::get_name | ( | ) | const [virtual] |
Return a string that uniquely identifies a finite element. This class returns FE_RaviartThomasNodal<dim>(degree)
, with dim
and degree
replaced by appropriate values.
Implements FiniteElement< dim, dim >.
virtual FiniteElement<dim>* FE_RaviartThomasNodal< dim >::clone | ( | ) | const [virtual] |
A sort of virtual copy constructor. Some places in the library, for example the constructors of FESystem as well as the hp::FECollection class, need to make copied of finite elements without knowing their exact type. They do so through this function.
Implements FiniteElement< dim, dim >.
virtual void FE_RaviartThomasNodal< dim >::interpolate | ( | std::vector< double > & | local_dofs, | |
const std::vector< double > & | values | |||
) | const [virtual] |
Interpolate a set of scalar values, computed in the generalized support points.
Reimplemented from FiniteElement< dim, dim >.
virtual void FE_RaviartThomasNodal< dim >::interpolate | ( | std::vector< double > & | local_dofs, | |
const std::vector< Vector< double > > & | values, | |||
unsigned int | offset = 0 | |||
) | const [virtual] |
Interpolate a set of vector values, computed in the generalized support points.
Since a finite element often only interpolates part of a vector, offset
is used to determine the first component of the vector to be interpolated. Maybe consider changing your data structures to use the next function.
Reimplemented from FiniteElement< dim, dim >.
virtual void FE_RaviartThomasNodal< dim >::interpolate | ( | std::vector< double > & | local_dofs, | |
const VectorSlice< const std::vector< std::vector< double > > > & | values | |||
) | const [virtual] |
Interpolate a set of vector values, computed in the generalized support points.
Reimplemented from FiniteElement< dim, dim >.
virtual void FE_RaviartThomasNodal< dim >::get_face_interpolation_matrix | ( | const FiniteElement< dim > & | source, | |
FullMatrix< double > & | matrix | |||
) | const [virtual] |
virtual void FE_RaviartThomasNodal< dim >::get_subface_interpolation_matrix | ( | const FiniteElement< dim > & | source, | |
const unsigned int | subface, | |||
FullMatrix< double > & | matrix | |||
) | const [virtual] |
virtual bool FE_RaviartThomasNodal< dim >::hp_constraints_are_implemented | ( | ) | const [virtual] |
Return whether this element implements its hanging node constraints in the new way, which has to be used to make elements "hp compatible". That means, the element properly implements the get_face_interpolation_matrix and get_subface_interpolation_matrix methods. Therefore the return value also indicates whether a call to the get_face_interpolation_matrix() method and the get_subface_interpolation_matrix() method will generate an error or not.
Currently the main purpose of this function is to allow the make_hanging_node_constraints method to decide whether the new procedures, which are supposed to work in the hp framework can be used, or if the old well verified but not hp capable functions should be used. Once the transition to the new scheme for computing the interface constraints is complete, this function will be superfluous and will probably go away.
Derived classes should implement this function accordingly. The default assumption is that a finite element does not provide hp capable face interpolation, and the default implementation therefore returns false
.
Reimplemented from FiniteElement< dim, dim >.
virtual std::vector<std::pair<unsigned int, unsigned int> > FE_RaviartThomasNodal< dim >::hp_vertex_dof_identities | ( | const FiniteElement< dim > & | fe_other | ) | const [virtual] |
virtual std::vector<std::pair<unsigned int, unsigned int> > FE_RaviartThomasNodal< dim >::hp_line_dof_identities | ( | const FiniteElement< dim > & | fe_other | ) | const [virtual] |
virtual std::vector<std::pair<unsigned int, unsigned int> > FE_RaviartThomasNodal< dim >::hp_quad_dof_identities | ( | const FiniteElement< dim > & | fe_other | ) | const [virtual] |
virtual FiniteElementDomination::Domination FE_RaviartThomasNodal< dim >::compare_for_face_domination | ( | const FiniteElement< dim > & | fe_other | ) | const [virtual] |
static std::vector<unsigned int> FE_RaviartThomasNodal< dim >::get_dpo_vector | ( | const unsigned int | degree | ) | [static, private] |
Only for internal use. Its full name is get_dofs_per_object_vector
function and it creates the dofs_per_object
vector that is needed within the constructor to be passed to the constructor of FiniteElementData
.
static std::vector<bool> FE_RaviartThomasNodal< dim >::get_ria_vector | ( | const unsigned int | degree | ) | [static, private] |
Compute the vector used for the restriction_is_additive
field passed to the base class's constructor.
virtual bool FE_RaviartThomasNodal< dim >::has_support_on_face | ( | const unsigned int | shape_index, | |
const unsigned int | face_index | |||
) | const [private, virtual] |
Check whether a shape function may be non-zero on a face.
Right now, this is only implemented for RT0 in 1D. Otherwise, returns always true
.
Reimplemented from FiniteElement< dim, dim >.
void FE_RaviartThomasNodal< dim >::initialize_support_points | ( | const unsigned int | rt_degree | ) | [private] |
Initialize the FiniteElement<dim>::generalized_support_points and FiniteElement<dim>::generalized_face_support_points fields. Called from the constructor.