|
| FE_Q_DG0 (const unsigned int p) |
|
| FE_Q_DG0 (const Quadrature< 1 > &points) |
|
virtual std::string | get_name () const override |
|
virtual void | convert_generalized_support_point_values_to_dof_values (const std::vector< Vector< double > > &support_point_values, std::vector< double > &nodal_values) const override |
|
virtual void | get_interpolation_matrix (const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const override |
|
virtual bool | has_support_on_face (const unsigned int shape_index, const unsigned int face_index) const override |
|
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > | get_constant_modes () const override |
|
virtual std::unique_ptr< FiniteElement< dim, spacedim > > | clone () const override |
|
virtual FiniteElementDomination::Domination | compare_for_domination (const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final |
|
virtual void | get_face_interpolation_matrix (const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override |
|
virtual void | get_subface_interpolation_matrix (const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix, const unsigned int face_no=0) const override |
|
virtual const FullMatrix< double > & | get_restriction_matrix (const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const override |
|
virtual const FullMatrix< double > & | get_prolongation_matrix (const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const override |
|
virtual unsigned int | face_to_cell_index (const unsigned int face_dof_index, const unsigned int face, const types::geometric_orientation combined_orientation=numbers::default_geometric_orientation) const override |
|
unsigned int | get_degree () const |
|
virtual UpdateFlags | requires_update_flags (const UpdateFlags update_flags) const override |
|
const ScalarPolynomialsBase< dim > & | get_poly_space () const |
|
std::vector< unsigned int > | get_poly_space_numbering () const |
|
std::vector< unsigned int > | get_poly_space_numbering_inverse () const |
|
virtual double | shape_value (const unsigned int i, const Point< dim > &p) const override |
|
virtual double | shape_value_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const override |
|
virtual Tensor< 1, dim > | shape_grad (const unsigned int i, const Point< dim > &p) const override |
|
virtual Tensor< 1, dim > | shape_grad_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const override |
|
virtual Tensor< 2, dim > | shape_grad_grad (const unsigned int i, const Point< dim > &p) const override |
|
virtual Tensor< 2, dim > | shape_grad_grad_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const override |
|
virtual Tensor< 3, dim > | shape_3rd_derivative (const unsigned int i, const Point< dim > &p) const override |
|
virtual Tensor< 3, dim > | shape_3rd_derivative_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const override |
|
virtual Tensor< 4, dim > | shape_4th_derivative (const unsigned int i, const Point< dim > &p) const override |
|
virtual Tensor< 4, dim > | shape_4th_derivative_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const override |
|
virtual std::size_t | memory_consumption () const override |
|
std::pair< std::unique_ptr< FiniteElement< dim, spacedim > >, unsigned int > | operator^ (const unsigned int multiplicity) const |
|
virtual bool | operator== (const FiniteElement< dim, spacedim > &fe) const |
|
bool | operator== (const FiniteElementData &) const |
|
bool | operator!= (const FiniteElement< dim, spacedim > &) const |
|
template<class Archive > |
void | serialize (Archive &ar, const unsigned int version) |
|
ReferenceCell | reference_cell () const |
|
unsigned int | n_unique_2d_subobjects () const |
|
unsigned int | n_unique_faces () const |
|
unsigned int | n_dofs_per_vertex () const |
|
unsigned int | n_dofs_per_line () const |
|
unsigned int | n_dofs_per_quad (unsigned int face_no=0) const |
|
unsigned int | max_dofs_per_quad () const |
|
unsigned int | n_dofs_per_hex () const |
|
unsigned int | n_dofs_per_face (unsigned int face_no=0, unsigned int child=0) const |
|
unsigned int | max_dofs_per_face () const |
|
unsigned int | n_dofs_per_cell () const |
|
template<int structdim> |
unsigned int | n_dofs_per_object (const unsigned int i=0) const |
|
unsigned int | n_components () const |
|
unsigned int | n_blocks () const |
|
const BlockIndices & | block_indices () const |
|
unsigned int | tensor_degree () const |
|
bool | conforms (const Conformity) const |
|
unsigned int | get_first_line_index () const |
|
unsigned int | get_first_quad_index (const unsigned int quad_no=0) const |
|
unsigned int | get_first_hex_index () const |
|
unsigned int | get_first_face_line_index (const unsigned int face_no=0) const |
|
unsigned int | get_first_face_quad_index (const unsigned int face_no=0) const |
|
|
virtual bool | hp_constraints_are_implemented () const override |
|
virtual std::vector< std::pair< unsigned int, unsigned int > > | hp_vertex_dof_identities (const FiniteElement< dim, spacedim > &fe_other) const override |
|
virtual std::vector< std::pair< unsigned int, unsigned int > > | hp_line_dof_identities (const FiniteElement< dim, spacedim > &fe_other) const override |
|
virtual std::vector< std::pair< unsigned int, unsigned int > > | hp_quad_dof_identities (const FiniteElement< dim, spacedim > &fe_other, const unsigned int face_no=0) const override |
|
|
bool | prolongation_is_implemented () const |
|
bool | isotropic_prolongation_is_implemented () const |
|
bool | restriction_is_implemented () const |
|
bool | isotropic_restriction_is_implemented () const |
|
bool | restriction_is_additive (const unsigned int index) const |
|
const FullMatrix< double > & | constraints (const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const |
|
bool | constraints_are_implemented (const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const |
|
|
std::pair< unsigned int, unsigned int > | system_to_component_index (const unsigned int index) const |
|
unsigned int | component_to_system_index (const unsigned int component, const unsigned int index) const |
|
std::pair< unsigned int, unsigned int > | face_system_to_component_index (const unsigned int index, const unsigned int face_no=0) const |
|
unsigned int | adjust_quad_dof_index_for_face_orientation (const unsigned int index, const unsigned int face_no, const types::geometric_orientation combined_orientation) const |
|
unsigned int | adjust_line_dof_index_for_line_orientation (const unsigned int index, const types::geometric_orientation combined_orientation) const |
|
const ComponentMask & | get_nonzero_components (const unsigned int i) const |
|
unsigned int | n_nonzero_components (const unsigned int i) const |
|
bool | is_primitive () const |
|
bool | is_primitive (const unsigned int i) const |
|
virtual const Table< 2, bool > & | get_local_dof_sparsity_pattern () const |
|
unsigned int | n_base_elements () const |
|
virtual const FiniteElement< dim, spacedim > & | base_element (const unsigned int index) const |
|
unsigned int | element_multiplicity (const unsigned int index) const |
|
const FiniteElement< dim, spacedim > & | get_sub_fe (const ComponentMask &mask) const |
|
virtual const FiniteElement< dim, spacedim > & | get_sub_fe (const unsigned int first_component, const unsigned int n_selected_components) const |
|
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > | system_to_base_index (const unsigned int index) const |
|
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > | face_system_to_base_index (const unsigned int index, const unsigned int face_no=0) const |
|
types::global_dof_index | first_block_of_base (const unsigned int b) const |
|
std::pair< unsigned int, unsigned int > | component_to_base_index (const unsigned int component) const |
|
std::pair< unsigned int, unsigned int > | block_to_base_index (const unsigned int block) const |
|
std::pair< unsigned int, types::global_dof_index > | system_to_block_index (const unsigned int component) const |
|
unsigned int | component_to_block_index (const unsigned int component) const |
|
bool | shape_function_belongs_to (const unsigned int shape_function, const FEValuesExtractors::Scalar &component) const |
|
bool | shape_function_belongs_to (const unsigned int shape_function, const FEValuesExtractors::Vector &components) const |
|
bool | shape_function_belongs_to (const unsigned int shape_function, const FEValuesExtractors::SymmetricTensor< 2 > &components) const |
|
bool | shape_function_belongs_to (const unsigned int shape_function, const FEValuesExtractors::Tensor< 2 > &components) const |
|
|
ComponentMask | component_mask (const FEValuesExtractors::Scalar &scalar) const |
|
ComponentMask | component_mask (const FEValuesExtractors::Vector &vector) const |
|
ComponentMask | component_mask (const FEValuesExtractors::SymmetricTensor< 2 > &sym_tensor) const |
|
ComponentMask | component_mask (const BlockMask &block_mask) const |
|
BlockMask | block_mask (const FEValuesExtractors::Scalar &scalar) const |
|
BlockMask | block_mask (const FEValuesExtractors::Vector &vector) const |
|
BlockMask | block_mask (const FEValuesExtractors::SymmetricTensor< 2 > &sym_tensor) const |
|
BlockMask | block_mask (const ComponentMask &component_mask) const |
|
|
const std::vector< Point< dim > > & | get_unit_support_points () const |
|
bool | has_support_points () const |
|
virtual Point< dim > | unit_support_point (const unsigned int index) const |
|
const std::vector< Point< dim - 1 > > & | get_unit_face_support_points (const unsigned int face_no=0) const |
|
bool | has_face_support_points (const unsigned int face_no=0) const |
|
virtual Point< dim - 1 > | unit_face_support_point (const unsigned int index, const unsigned int face_no=0) const |
|
const std::vector< Point< dim > > & | get_generalized_support_points () const |
|
bool | has_generalized_support_points () const |
|
GeometryPrimitive | get_associated_geometry_primitive (const unsigned int cell_dof_index) const |
|
|
unsigned int | n_subscriptions () const |
|
template<typename StreamType > |
void | list_subscribers (StreamType &stream) const |
|
void | list_subscribers () const |
|
|
void | initialize (const std::vector< Point< 1 > > &support_points_1d) |
|
void | initialize_constraints (const std::vector< Point< 1 > > &points) |
|
void | initialize_unit_support_points (const std::vector< Point< 1 > > &points) |
|
void | initialize_unit_face_support_points (const std::vector< Point< 1 > > &points) |
|
void | initialize_dof_index_permutations () |
|
virtual std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > | get_data (const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override |
|
virtual void | fill_fe_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override |
|
virtual void | fill_fe_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const =0 |
|
virtual void | fill_fe_face_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override |
|
virtual void | fill_fe_face_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const |
|
virtual void | fill_fe_face_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const |
|
virtual void | fill_fe_subface_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override |
|
virtual void | fill_fe_subface_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const =0 |
|
void | correct_hessians (internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const unsigned int n_q_points) const |
|
void | correct_third_derivatives (internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const unsigned int n_q_points) const |
|
void | reinit_restriction_and_prolongation_matrices (const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false) |
|
TableIndices< 2 > | interface_constraints_size () const |
|
virtual std::unique_ptr< InternalDataBase > | get_face_data (const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const hp::QCollection< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const |
|
virtual std::unique_ptr< InternalDataBase > | get_face_data (const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim - 1 > &quadrature, internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const |
|
virtual std::unique_ptr< InternalDataBase > | get_subface_data (const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const |
|
|
const std::unique_ptr< ScalarPolynomialsBase< dim > > | poly_space |
|
std::vector< std::vector< FullMatrix< double > > > | restriction |
|
std::vector< std::vector< FullMatrix< double > > > | prolongation |
|
FullMatrix< double > | interface_constraints |
|
std::vector< Point< dim > > | unit_support_points |
|
std::vector< std::vector< Point< dim - 1 > > > | unit_face_support_points |
|
std::vector< Point< dim > > | generalized_support_points |
|
std::vector< std::vector< Point< dim - 1 > > > | generalized_face_support_points |
|
std::vector< Table< 2, int > > | adjust_quad_dof_index_for_face_orientation_table |
|
std::vector< int > | adjust_line_dof_index_for_line_orientation_table |
|
std::vector< std::pair< unsigned int, unsigned int > > | system_to_component_table |
|
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > | face_system_to_component_table |
|
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > | system_to_base_table |
|
std::vector< std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > > | face_system_to_base_table |
|
BlockIndices | base_to_block_indices |
|
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > | component_to_base_table |
|
const std::vector< bool > | restriction_is_additive_flags |
|
const std::vector< ComponentMask > | nonzero_components |
|
const std::vector< unsigned int > | n_nonzero_components_table |
|
const bool | cached_primitivity |
|
Table< 2, bool > | local_dof_sparsity_pattern |
|
template<
int dim,
int spacedim = dim>
class FE_Q_DG0< dim, spacedim >
Implementation of a scalar Lagrange finite element Qp+DG0
that yields the finite element space of continuous, piecewise polynomials of degree p
in each coordinate direction plus the space of locally constant functions. This class is realized using tensor product polynomials based on equidistant or given support points.
The standard constructor of this class takes the degree p
of this finite element. Alternatively, it can take a quadrature formula points
defining the support points of the Lagrange interpolation in one coordinate direction.
For more information about the spacedim
template parameter check the documentation of FiniteElement or the one of Triangulation.
For more information regarding this element see [37] .
Implementation
The constructor creates a TensorProductPolynomials object that includes the tensor product of LagrangeEquidistant
polynomials of degree p
plus the locally constant function. This TensorProductPolynomialsConst
object provides all values and derivatives of the shape functions. In case a quadrature rule is given, the constructor creates a TensorProductPolynomialsConst object that includes the tensor product of Lagrange
polynomials with the support points from points
and a locally constant function.
Furthermore the constructor fills the interface_constrains
, the prolongation
(embedding) and the restriction
matrices.
Numbering of the degrees of freedom (DoFs)
The original ordering of the shape functions represented by the TensorProductPolynomialsConst is a tensor product numbering. However, the shape functions on a cell are renumbered beginning with the shape functions whose support points are at the vertices, then on the line, on the quads, and finally (for 3d) on the hexes. Finally there is a support point for the discontinuous shape function in the middle of the cell. To be explicit, these numberings are listed in the following:
Q1 elements
-
1d case:
* 0---2---1
*
-
2d case:
* 2-------3
* | |
* | 5 |
* | |
* 0-------1
*
-
3d case:
* 6-------7 6-------7
* /| | / /|
* / | | / / |
* / | | / / |
* 4 | 8 | 4-------5 |
* | 2-------3 | | 3
* | / / | | /
* | / / | | /
* |/ / | |/
* 0-------1 0-------1
*
The respective coordinate values of the support points of the degrees of freedom are as follows:
-
Index 0:
[ 0, 0, 0]
;
-
Index 1:
[ 1, 0, 0]
;
-
Index 2:
[ 0, 1, 0]
;
-
Index 3:
[ 1, 1, 0]
;
-
Index 4:
[ 0, 0, 1]
;
-
Index 5:
[ 1, 0, 1]
;
-
Index 6:
[ 0, 1, 1]
;
-
Index 7:
[ 1, 1, 1]
;
-
Index 8:
[1/2, 1/2, 1/2]
;
Q2 elements
-
1d case:
* 0---2---1
*
Index 3 has the same coordinates as index 2
-
2d case:
* 2---7---3
* | |
* 4 8 5
* | |
* 0---6---1
*
Index 9 has the same coordinates as index 2
-
3d case:
* 6--15---7 6--15---7
* /| | / /|
* 12 | 19 12 1319
* / 18 | / / |
* 4 | | 4---14--5 |
* | 2---11--3 | | 3
* | / / | 17 /
* 16 8 9 16 | 9
* |/ / | |/
* 0---10--1 0---8---1
*
* *-------* *-------*
* /| | / /|
* / | 23 | / 25 / |
* / | | / / |
* * | | *-------* |
* |20 *-------* | |21 *
* | / / | 22 | /
* | / 24 / | | /
* |/ / | |/
* *-------* *-------*
*
The center vertices have number 26 and 27.
The respective coordinate values of the support points of the degrees of freedom are as follows:
-
Index 0:
[0, 0, 0]
;
-
Index 1:
[1, 0, 0]
;
-
Index 2:
[0, 1, 0]
;
-
Index 3:
[1, 1, 0]
;
-
Index 4:
[0, 0, 1]
;
-
Index 5:
[1, 0, 1]
;
-
Index 6:
[0, 1, 1]
;
-
Index 7:
[1, 1, 1]
;
-
Index 8:
[0, 1/2, 0]
;
-
Index 9:
[1, 1/2, 0]
;
-
Index 10:
[1/2, 0, 0]
;
-
Index 11:
[1/2, 1, 0]
;
-
Index 12:
[0, 1/2, 1]
;
-
Index 13:
[1, 1/2, 1]
;
-
Index 14:
[1/2, 0, 1]
;
-
Index 15:
[1/2, 1, 1]
;
-
Index 16:
[0, 0, 1/2]
;
-
Index 17:
[1, 0, 1/2]
;
-
Index 18:
[0, 1, 1/2]
;
-
Index 19:
[1, 1, 1/2]
;
-
Index 20:
[0, 1/2, 1/2]
;
-
Index 21:
[1, 1/2, 1/2]
;
-
Index 22:
[1/2, 0, 1/2]
;
-
Index 23:
[1/2, 1, 1/2]
;
-
Index 24:
[1/2, 1/2, 0]
;
-
Index 25:
[1/2, 1/2, 1]
;
-
Index 26:
[1/2, 1/2, 1/2]
;
-
Index 27:
[1/2, 1/2, 1/2]
;
Q3 elements
Q4 elements
Definition at line 237 of file fe_q_dg0.h.
template<
int dim,
int spacedim>
void FE_Q_DG0< dim, spacedim >::convert_generalized_support_point_values_to_dof_values |
( |
const std::vector< Vector< double > > & |
support_point_values, |
|
|
std::vector< double > & |
nodal_values |
|
) |
| const |
|
overridevirtual |
Given the values of a function \(f(\mathbf x)\) at the (generalized) support points of the reference cell, this function then computes what the nodal values of the element are, i.e., \(\Psi_i[f]\), where \(\Psi_i\) are the node functionals of the element (see also Node values or node functionals). The values \(\Psi_i[f]\) are then the expansion coefficients for the shape functions of the finite element function that interpolates the given function \(f(x)\), i.e., \( f_h(\mathbf x) = \sum_i \Psi_i[f] \varphi_i(\mathbf x)
\) is the finite element interpolant of \(f\) with the current element. The operation described here is used, for example, in the FETools::compute_node_matrix() function.
In more detail, let us assume that the generalized support points (see this glossary entry ) of the current element are \(\hat{\mathbf x}_i\) and that the node functionals associated with the current element are \(\Psi_i[\cdot]\). Then, the fact that the element is based on generalized support points, implies that if we apply \(\Psi_i\) to a (possibly vector-valued) finite element function \(\varphi\), the result must have the form \(\Psi_i[\varphi] = f_i(\varphi(\hat{\mathbf x}_i))\) – in other words, the value of the node functional \(\Psi_i\) applied to \(\varphi\) only depends on the values of \(\varphi\) at \(\hat{\mathbf x}_i\) and not on values anywhere else, or integrals of \(\varphi\), or any other kind of information.
The exact form of \(f_i\) depends on the element. For example, for scalar Lagrange elements, we have that in fact \(\Psi_i[\varphi] = \varphi(\hat{\mathbf x}_i)\). If you combine multiple scalar Lagrange elements via an FESystem object, then \(\Psi_i[\varphi] = \varphi(\hat{\mathbf x}_i)_{c(i)}\) where \(c(i)\) is the result of the FiniteElement::system_to_component_index() function's return value's first component. In these two cases, \(f_i\) is therefore simply the identity (in the scalar case) or a function that selects a particular vector component of its argument. On the other hand, for Raviart-Thomas elements, one would have that \(f_i(\mathbf y) = \mathbf y \cdot \mathbf n_i\) where \(\mathbf n_i\) is the normal vector of the face at which the shape function is defined.
Given all of this, what this function does is the following: If you input a list of values of a function \(\varphi\) at all generalized support points (where each value is in fact a vector of values with as many components as the element has), then this function returns a vector of values obtained by applying the node functionals to these values. In other words, if you pass in \(\{\varphi(\hat{\mathbf x}_i)\}_{i=0}^{N-1}\) then you will get out a vector \(\{\Psi[\varphi]\}_{i=0}^{N-1}\) where \(N\) equals dofs_per_cell
.
- Parameters
-
[in] | support_point_values | An array of size dofs_per_cell (which equals the number of points the get_generalized_support_points() function will return) where each element is a vector with as many entries as the element has vector components. This array should contain the values of a function at the generalized support points of the current element. |
[out] | nodal_values | An array of size dofs_per_cell that contains the node functionals of the element applied to the given function. |
- Note
- It is safe to call this function for (transformed) values on the real cell only for elements with trivial MappingKind. For all other elements (for example for H(curl), or H(div) conforming elements) vector values have to be transformed to the reference cell first.
-
Given what the function is supposed to do, the function clearly can only work for elements that actually implement (generalized) support points. Elements that do not have generalized support points – e.g., elements whose nodal functionals evaluate integrals or moments of functions (such as FE_Q_Hierarchical) – can in general not make sense of the operation that is required for this function. They consequently may not implement it.
Reimplemented from FiniteElement< dim, spacedim >.
Definition at line 180 of file fe_q_dg0.cc.
template<
int dim,
int spacedim = dim>
Projection from a fine grid space onto a coarse grid space. Overrides the respective method in FiniteElement, implementing lazy evaluation (initialize when requested).
If this projection operator is associated with a matrix P
, then the restriction of this matrix P_i
to a single child cell is returned here.
The matrix P
is the concatenation or the sum of the cell matrices P_i
, depending on the restriction_is_additive_flags. This distinguishes interpolation (concatenation) and projection with respect to scalar products (summation).
Row and column indices are related to coarse grid and fine grid spaces, respectively, consistent with the definition of the associated operator.
If projection matrices are not implemented in the derived finite element class, this function aborts with ExcProjectionVoid. You can check whether this is the case by calling the restriction_is_implemented() or the isotropic_restriction_is_implemented() function.
Reimplemented from FiniteElement< dim, spacedim >.
Reimplemented in FE_Q_Bubbles< dim, spacedim >, and FE_Bernstein< dim, spacedim >.
template<
int dim,
int spacedim = dim>
Embedding matrix between grids. Overrides the respective method in FiniteElement, implementing lazy evaluation (initialize when queried).
The identity operator from a coarse grid space into a fine grid space is associated with a matrix P
. The restriction of this matrix P_i
to a single child cell is returned here.
The matrix P
is the concatenation, not the sum of the cell matrices P_i
. That is, if the same non-zero entry j,k
exists in two different child matrices P_i
, the value should be the same in both matrices and it is copied into the matrix P
only once.
Row and column indices are related to fine grid and coarse grid spaces, respectively, consistent with the definition of the associated operator.
These matrices are used by routines assembling the prolongation matrix for multi-level methods. Upon assembling the transfer matrix between cells using this matrix array, zero elements in the prolongation matrix are discarded and will not fill up the transfer matrix.
If projection matrices are not implemented in the derived finite element class, this function aborts with ExcEmbeddingVoid. You can check whether this is the case by calling the prolongation_is_implemented() or the isotropic_prolongation_is_implemented() function.
Reimplemented from FiniteElement< dim, spacedim >.
Reimplemented in FE_Q_Bubbles< dim, spacedim >, and FE_Bernstein< dim, spacedim >.
template<
int dim,
int spacedim = dim>
|
inlineoverrideprotectedvirtualinherited |
Create an internal data object and return a pointer to it of which the caller of this function then assumes ownership. This object will then be passed to the FiniteElement::fill_fe_values() every time the finite element shape functions and their derivatives are evaluated on a concrete cell. The object created here is therefore used by derived classes as a place for scratch objects that are used in evaluating shape functions, as well as to store information that can be pre-computed once and re-used on every cell (e.g., for evaluating the values and gradients of shape functions on the reference cell, for later re-use when transforming these values to a concrete cell).
This function is the first one called in the process of initializing a FEValues object for a given mapping and finite element object. The returned object will later be passed to FiniteElement::fill_fe_values() for a concrete cell, which will itself place its output into an object of type internal::FEValuesImplementation::FiniteElementRelatedData. Since there may be data that can already be computed in its final form on the reference cell, this function also receives a reference to the internal::FEValuesImplementation::FiniteElementRelatedData object as its last argument. This output argument is guaranteed to always be the same one when used with the InternalDataBase object returned by this function. In other words, the subdivision of scratch data and final data in the returned object and the output_data
object is as follows: If data can be pre-computed on the reference cell in the exact form in which it will later be needed on a concrete cell, then this function should already emplace it in the output_data
object. An example are the values of shape functions at quadrature points for the usual Lagrange elements which on a concrete cell are identical to the ones on the reference cell. On the other hand, if some data can be pre-computed to make computations on a concrete cell cheaper, then it should be put into the returned object for later re-use in a derive class's implementation of FiniteElement::fill_fe_values(). An example are the gradients of shape functions on the reference cell for Lagrange elements: to compute the gradients of the shape functions on a concrete cell, one has to multiply the gradients on the reference cell by the inverse of the Jacobian of the mapping; consequently, we cannot already compute the gradients on a concrete cell at the time the current function is called, but we can at least pre-compute the gradients on the reference cell, and store it in the object returned.
An extensive discussion of the interaction between this function and FEValues can be found in the How Mapping, FiniteElement, and FEValues work together documentation topic. See also the documentation of the InternalDataBase class.
- Parameters
-
[in] | update_flags | A set of UpdateFlags values that describe what kind of information the FEValues object requests the finite element to compute. This set of flags may also include information that the finite element can not compute, e.g., flags that pertain to data produced by the mapping. An implementation of this function needs to set up all data fields in the returned object that are necessary to produce the finite-element related data specified by these flags, and may already pre-compute part of this information as discussed above. Elements may want to store these update flags (or a subset of these flags) in InternalDataBase::update_each so they know at the time when FiniteElement::fill_fe_values() is called what they are supposed to compute |
[in] | mapping | A reference to the mapping used for computing values and derivatives of shape functions. |
[in] | quadrature | A reference to the object that describes where the shape functions should be evaluated. |
[out] | output_data | A reference to the object that FEValues will use in conjunction with the object returned here and where an implementation of FiniteElement::fill_fe_values() will place the requested information. This allows the current function to already pre-compute pieces of information that can be computed on the reference cell, as discussed above. FEValues guarantees that this output object and the object returned by the current function will always be used together. |
- Returns
- A pointer to an object of a type derived from InternalDataBase and that derived classes can use to store scratch data that can be pre-computed, or for scratch arrays that then only need to be allocated once. The calling site assumes ownership of this object and will delete it when it is no longer necessary.
Implements FiniteElement< dim, spacedim >.
Reimplemented in FE_Hermite< dim, spacedim >.
Definition at line 262 of file fe_poly.h.
template<
int dim,
int spacedim = dim>
virtual void FiniteElement< dim, spacedim >::fill_fe_values |
( |
const typename Triangulation< dim, spacedim >::cell_iterator & |
cell, |
|
|
const CellSimilarity::Similarity |
cell_similarity, |
|
|
const Quadrature< dim > & |
quadrature, |
|
|
const Mapping< dim, spacedim > & |
mapping, |
|
|
const typename Mapping< dim, spacedim >::InternalDataBase & |
mapping_internal, |
|
|
const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > & |
mapping_data, |
|
|
const InternalDataBase & |
fe_internal, |
|
|
::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > & |
output_data |
|
) |
| const |
|
protectedpure virtualinherited |
Compute information about the shape functions on the cell denoted by the first argument. Derived classes will have to implement this function based on the kind of element they represent. It is called by FEValues::reinit().
Conceptually, this function evaluates shape functions and their derivatives at the quadrature points represented by the mapped locations of those described by the quadrature argument to this function. In many cases, computing derivatives of shape functions (and in some cases also computing values of shape functions) requires making use of the mapping from the reference to the real cell; this information can either be taken from the mapping_data
object that has been filled for the current cell before this function is called, or by calling the member functions of a Mapping object with the mapping_internal
object that also corresponds to the current cell.
The information computed by this function is used to fill the various member variables of the output argument of this function. Which of the member variables of that structure should be filled is determined by the update flags stored in the FiniteElement::InternalDataBase::update_each field of the object passed to this function. These flags are typically set by FiniteElement::get_data(), FiniteElement::get_face_date() and FiniteElement::get_subface_data() (or, more specifically, implementations of these functions in derived classes).
An extensive discussion of the interaction between this function and FEValues can be found in the How Mapping, FiniteElement, and FEValues work together documentation topic.
- Parameters
-
[in] | cell | The cell of the triangulation for which this function is to compute a mapping from the reference cell to. |
[in] | cell_similarity | Whether or not the cell given as first argument is simply a translation, rotation, etc of the cell for which this function was called the most recent time. This information is computed simply by matching the vertices (as stored by the Triangulation) between the previous and the current cell. The value passed here may be modified by implementations of this function and should then be returned (see the discussion of the return value of this function). |
[in] | quadrature | A reference to the quadrature formula in use for the current evaluation. This quadrature object is the same as the one used when creating the internal_data object. The current object is then responsible for evaluating shape functions at the mapped locations of the quadrature points represented by this object. |
[in] | mapping | A reference to the mapping object used to map from the reference cell to the current cell. This object was used to compute the information in the mapping_data object before the current function was called. It is also the mapping object that created the mapping_internal object via Mapping::get_data(). You will need the reference to this mapping object most often to call Mapping::transform() to transform gradients and higher derivatives from the reference to the current cell. |
[in] | mapping_internal | An object specific to the mapping object. What the mapping chooses to store in there is of no relevance to the current function, but you may have to pass a reference to this object to certain functions of the Mapping class (e.g., Mapping::transform()) if you need to call them from the current function. |
[in] | mapping_data | The output object into which the Mapping::fill_fe_values() function wrote the mapping information corresponding to the current cell. This includes, for example, Jacobians of the mapping that may be of relevance to the current function, as well as other information that FEValues::reinit() requested from the mapping. |
[in] | fe_internal | A reference to an object previously created by get_data() and that may be used to store information the mapping can compute once on the reference cell. See the documentation of the FiniteElement::InternalDataBase class for an extensive description of the purpose of these objects. |
[out] | output_data | A reference to an object whose member variables should be computed. Not all of the members of this argument need to be filled; which ones need to be filled is determined by the update flags stored inside the fe_internal object. |
- Note
- FEValues ensures that this function is always called with the same pair of
fe_internal
and output_data
objects. In other words, if an implementation of this function knows that it has written a piece of data into the output argument in a previous call, then there is no need to copy it there again in a later call if the implementation knows that this is the same value.
template<
int dim,
int spacedim = dim>
Compute vector component and index of this shape function within the shape functions corresponding to this component from the index of a shape function within this finite element.
If the element is scalar, then the component is always zero, and the index within this component is equal to the overall index.
If the shape function referenced has more than one non-zero component, then it cannot be associated with one vector component, and an exception of type ExcShapeFunctionNotPrimitive will be raised.
Note that if the element is composed of other (base) elements, and a base element has more than one component but all its shape functions are primitive (i.e. are non-zero in only one component), then this mapping contains valid information. However, the index of a shape function of this element within one component (i.e. the second number of the respective entry of this array) does not indicate the index of the respective shape function within the base element (since that has more than one vector-component). For this information, refer to the system_to_base_table field and the system_to_base_index() function.
See the class description above for an example of how this function is typically used.
The use of this function is explained extensively in the step-8 and step-20 tutorial programs as well as in the Handling vector valued problems topic.
template<
int dim,
int spacedim = dim>
Return whether a finite element has defined support points. If the result is true, then a call to the get_unit_support_points() yields an array with dofs_per_cell
entries. Note that the function's name is poorly chosen: The function does not return whether an element has support points, but whether its implementation is able to provide a list of support points via get_unit_support_points().
The result may be false if a finite element defines itself not by interpolating shape functions, but by other means. A typical example are discontinuous \(P\)-type elements on quadrilaterals (rather than the common \(Q\)-type elements on quadrilaterals). Elements will generally only return true
if they construct their shape functions by the requirement that they be nonzero at a certain point and zero at all the points associated with the other shape functions. In other words, if the "node functionals" that form the dual space of the finite element space and are used to define the shape functions are point evaluations.
In composed elements (i.e. for the FESystem class), the result will be true if all the base elements have defined support points. FE_Nothing is a special case since it does not have support point (it has no shape functions, after all), and so returns an empty array from its get_unit_support_points() function. Nonetheless, because that array has the right size, it reports true
in the current function.
template<
int dim,
int spacedim = dim>
Return the support points of the trial functions on the unit face, if the derived finite element defines some. Finite elements that allow some kind of interpolation operation usually have support points. On the other hand, elements that define their degrees of freedom by, for example, moments on faces, or as derivatives, don't have support points. In that case, the returned field is empty
Note that elements that have support points need not necessarily have some on the faces, even if the interpolation points are located physically on a face. For example, the discontinuous elements have interpolation points on the vertices, and for higher degree elements also on the faces, but they are not defined to be on faces since in that case degrees of freedom from both sides of a face (or from all adjacent elements to a vertex) would be identified with each other, which is not what we would like to have). Logically, these degrees of freedom are therefore defined to belong to the cell, rather than the face or vertex. In that case, the returned element would therefore have length zero.
If the finite element defines support points, then their number equals the number of degrees of freedom on the face (dofs_per_face). The order of points in the array matches that returned by the cell->face(face)->get_dof_indices
function.
See the class documentation for details on support points.
template<
int dim,
int spacedim = dim>
For a given degree of freedom, return whether it is logically associated with a vertex, line, quad or hex.
For instance, for continuous finite elements this coincides with the lowest dimensional object the support point of the degree of freedom lies on. To give an example, for \(Q_1\) elements in 3d, every degree of freedom is defined by a shape function that we get by interpolating using support points that lie on the vertices of the cell. The support of these points of course extends to all edges connected to this vertex, as well as the adjacent faces and the cell interior, but we say that logically the degree of freedom is associated with the vertex as this is the lowest-dimensional object it is associated with. Likewise, for \(Q_2\) elements in 3d, the degrees of freedom with support points at edge midpoints would yield a value of GeometryPrimitive::line from this function, whereas those on the centers of faces in 3d would return GeometryPrimitive::quad.
To make this more formal, the kind of object returned by this function represents the object so that the support of the shape function corresponding to the degree of freedom, (i.e., that part of the domain where the function "lives") is the union of all of the cells sharing this object. To return to the example above, for \(Q_2\) in 3d, the shape function with support point at an edge midpoint has support on all cells that share the edge and not only the cells that share the adjacent faces, and consequently the function will return GeometryPrimitive::line.
On the other hand, for discontinuous elements of type \(DGQ_2\), a degree of freedom associated with an interpolation polynomial that has its support point physically located at a line bounding a cell, but is nonzero only on one cell. Consequently, it is logically associated with the interior of that cell (i.e., with a GeometryPrimitive::quad in 2d and a GeometryPrimitive::hex in 3d).
- Parameters
-
[in] | cell_dof_index | The index of a shape function or degree of freedom. This index must be in the range [0,dofs_per_cell) . |
- Note
- The integer value of the object returned by this function equals the dimensionality of the object it describes, and can consequently be used in generic programming paradigms. For example, if a degree of freedom is associated with a vertex, then this function returns GeometryPrimitive::vertex, which has a numeric value of zero (the dimensionality of a vertex).
template<
int dim,
int spacedim = dim>
|
protectedvirtualinherited |
Like get_data(), but return an object that will later be used for evaluating shape function information at quadrature points on faces of cells. The object will then be used in calls to implementations of FiniteElement::fill_fe_face_values(). See the documentation of get_data() for more information.
The default implementation of this function converts the face quadrature into a cell quadrature with appropriate quadrature point locations, and with that calls the get_data() function above that has to be implemented in derived classes.
- Parameters
-
[in] | update_flags | A set of UpdateFlags values that describe what kind of information the FEValues object requests the finite element to compute. This set of flags may also include information that the finite element can not compute, e.g., flags that pertain to data produced by the mapping. An implementation of this function needs to set up all data fields in the returned object that are necessary to produce the finite-element related data specified by these flags, and may already pre-compute part of this information as discussed above. Elements may want to store these update flags (or a subset of these flags) in InternalDataBase::update_each so they know at the time when FiniteElement::fill_fe_face_values() is called what they are supposed to compute |
[in] | mapping | A reference to the mapping used for computing values and derivatives of shape functions. |
[in] | quadrature | A reference to the object that describes where the shape functions should be evaluated. |
[out] | output_data | A reference to the object that FEValues will use in conjunction with the object returned here and where an implementation of FiniteElement::fill_fe_face_values() will place the requested information. This allows the current function to already pre-compute pieces of information that can be computed on the reference cell, as discussed above. FEValues guarantees that this output object and the object returned by the current function will always be used together. |
- Returns
- A pointer to an object of a type derived from InternalDataBase and that derived classes can use to store scratch data that can be pre-computed, or for scratch arrays that then only need to be allocated once. The calling site assumes ownership of this object and will delete it when it is no longer necessary.
Reimplemented in FE_Enriched< dim, spacedim >, FE_PolyFace< PolynomialType, dim, spacedim >, FE_PolyFace< PolynomialSpace< dim - 1 >, dim, dim >, FE_PolyFace< TensorProductPolynomials< dim - 1 >, dim, dim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.
template<
int dim,
int spacedim = dim>
|
protectedvirtualinherited |
Like get_data(), but return an object that will later be used for evaluating shape function information at quadrature points on children of faces of cells. The object will then be used in calls to implementations of FiniteElement::fill_fe_subface_values(). See the documentation of get_data() for more information.
The default implementation of this function converts the face quadrature into a cell quadrature with appropriate quadrature point locations, and with that calls the get_data() function above that has to be implemented in derived classes.
- Parameters
-
[in] | update_flags | A set of UpdateFlags values that describe what kind of information the FEValues object requests the finite element to compute. This set of flags may also include information that the finite element can not compute, e.g., flags that pertain to data produced by the mapping. An implementation of this function needs to set up all data fields in the returned object that are necessary to produce the finite-element related data specified by these flags, and may already pre-compute part of this information as discussed above. Elements may want to store these update flags (or a subset of these flags) in InternalDataBase::update_each so they know at the time when FiniteElement::fill_fe_subface_values() is called what they are supposed to compute |
[in] | mapping | A reference to the mapping used for computing values and derivatives of shape functions. |
[in] | quadrature | A reference to the object that describes where the shape functions should be evaluated. |
[out] | output_data | A reference to the object that FEValues will use in conjunction with the object returned here and where an implementation of FiniteElement::fill_fe_subface_values() will place the requested information. This allows the current function to already pre-compute pieces of information that can be computed on the reference cell, as discussed above. FEValues guarantees that this output object and the object returned by the current function will always be used together. |
- Returns
- A pointer to an object of a type derived from InternalDataBase and that derived classes can use to store scratch data that can be pre-computed, or for scratch arrays that then only need to be allocated once. The calling site assumes ownership of this object and will delete it when it is no longer necessary.
Reimplemented in FE_Enriched< dim, spacedim >, FE_PolyFace< PolynomialType, dim, spacedim >, FE_PolyFace< PolynomialSpace< dim - 1 >, dim, dim >, FE_PolyFace< TensorProductPolynomials< dim - 1 >, dim, dim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.