|
| FE_FaceQ (const unsigned int p) |
|
virtual std::unique_ptr< FiniteElement< dim, spacedim > > | clone () const override |
|
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_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 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 |
|
unsigned int | get_degree () const |
|
virtual UpdateFlags | requires_update_flags (const UpdateFlags update_flags) 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 |
|
virtual std::size_t | memory_consumption () 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 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 |
|
virtual bool | hp_constraints_are_implemented () const override |
|
virtual FiniteElementDomination::Domination | compare_for_domination (const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final |
|
|
virtual double | shape_value (const unsigned int i, const Point< dim > &p) const |
|
virtual double | shape_value_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const |
|
virtual Tensor< 1, dim > | shape_grad (const unsigned int i, const Point< dim > &p) const |
|
virtual Tensor< 1, dim > | shape_grad_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const |
|
virtual Tensor< 2, dim > | shape_grad_grad (const unsigned int i, const Point< dim > &p) const |
|
virtual Tensor< 2, dim > | shape_grad_grad_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const |
|
virtual Tensor< 3, dim > | shape_3rd_derivative (const unsigned int i, const Point< dim > &p) const |
|
virtual Tensor< 3, dim > | shape_3rd_derivative_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const |
|
virtual Tensor< 4, dim > | shape_4th_derivative (const unsigned int i, const Point< dim > &p) const |
|
virtual Tensor< 4, dim > | shape_4th_derivative_component (const unsigned int i, const Point< dim > &p, const unsigned int component) const |
|
|
virtual const FullMatrix< double > & | get_restriction_matrix (const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const |
|
virtual const FullMatrix< double > & | get_prolongation_matrix (const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const |
|
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 |
|
virtual void | get_interpolation_matrix (const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) 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 |
|
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 |
|
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 |
|
|
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 |
|
std::unique_ptr< typename FiniteElement< dim, spacedim >::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 override |
|
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 |
|
std::unique_ptr< typename FiniteElement< dim, spacedim >::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 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 | reinit_restriction_and_prolongation_matrices (const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false) |
|
TableIndices< 2 > | interface_constraints_size () const |
|
|
PolynomialType | 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_FaceQ< dim, spacedim >
A finite element that is a tensor product polynomial on each face and undefined in the interior of the cells. The basis functions on the faces are Lagrange polynomials based on the support points of the (dim-1)-dimensional Gauss–Lobatto quadrature rule. For element degree one and two, the polynomials hence correspond to the usual Lagrange polynomials on equidistant points.
Although the name does not give it away, the element is discontinuous at locations where faces of cells meet. In particular, this finite element is the trace space of FE_RaviartThomas on the faces and serves in hybridized methods, e.g. in combination with the FE_DGQ element. Its use is demonstrated in the step-51 tutorial program.
- Note
- Since this element is defined only on faces, only FEFaceValues and FESubfaceValues will provide useful information. On the other hand, if you use this element with FEValues for cell integration, then the values and derivatives of shape functions will have invalid values and will not likely produce anything useful. In order to make the use of this element as part of an FESystem simpler, using a (cell) FEValues object will not fail outright, but those components of shape functions of the combined element that correspond to FE_FaceQ will have the invalid values mentioned above.
Definition at line 54 of file fe_face.h.
template<
typename PolynomialType ,
int dim = PolynomialType::dimension + 1,
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 >.
Definition at line 91 of file fe_poly_face.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>
Return the value of the ith
shape function at the point p
. p
is a point on the reference element. If the finite element is vector-valued, then return the value of the only non-zero component of the vector value of this shape function. If the shape function has more than one non-zero component (which we refer to with the term non-primitive), then derived classes implementing this function should throw an exception of type ExcShapeFunctionNotPrimitive. In that case, use the shape_value_component() function.
Implementations of this function should throw an exception of type ExcUnitShapeValuesDoNotExist if the shape functions of the FiniteElement under consideration depend on the shape of the cell in real space, i.e., if the shape functions are not defined by mapping from the reference cell. Some non-conforming elements are defined this way, as is the FE_DGPNonparametric class, to name just one example.
The default implementation of this virtual function does exactly this, i.e., it simply throws an exception of type ExcUnitShapeValuesDoNotExist.
Reimplemented in FE_DGPNonparametric< dim, spacedim >, FE_Enriched< dim, spacedim >, FE_NedelecSZ< dim, spacedim >, FE_Nothing< dim, spacedim >, FE_Nothing< dim, dim >, FE_Poly< dim, spacedim >, FE_Poly< dim, dim >, FE_Poly< dim, spacedim >, FE_PolyTensor< dim, spacedim >, FE_PolyTensor< dim, dim >, FE_PolyTensor< dim, spacedim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.
template<
int dim,
int spacedim = dim>
Return the gradient of the ith
shape function at the point p
. p
is a point on the reference element, and likewise the gradient is the gradient on the unit cell with respect to unit cell coordinates. If the finite element is vector-valued, then return the value of the only non-zero component of the vector value of this shape function. If the shape function has more than one non-zero component (which we refer to with the term non-primitive), then derived classes implementing this function should throw an exception of type ExcShapeFunctionNotPrimitive. In that case, use the shape_grad_component() function.
Implementations of this function should throw an exception of type ExcUnitShapeValuesDoNotExist if the shape functions of the FiniteElement under consideration depend on the shape of the cell in real space, i.e., if the shape functions are not defined by mapping from the reference cell. Some non-conforming elements are defined this way, as is the FE_DGPNonparametric class, to name just one example.
The default implementation of this virtual function does exactly this, i.e., it simply throws an exception of type ExcUnitShapeValuesDoNotExist.
Reimplemented in FE_DGPNonparametric< dim, spacedim >, FE_NedelecSZ< dim, spacedim >, FE_Poly< dim, spacedim >, FE_Poly< dim, dim >, FE_Poly< dim, spacedim >, FE_PolyTensor< dim, spacedim >, FE_PolyTensor< dim, dim >, FE_PolyTensor< dim, spacedim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.
template<
int dim,
int spacedim = dim>
Return the tensor of second derivatives of the ith
shape function at point p
on the unit cell. The derivatives are derivatives on the unit cell with respect to unit cell coordinates. If the finite element is vector-valued, then return the value of the only non-zero component of the vector value of this shape function. If the shape function has more than one non-zero component (which we refer to with the term non-primitive), then derived classes implementing this function should throw an exception of type ExcShapeFunctionNotPrimitive. In that case, use the shape_grad_grad_component() function.
Implementations of this function should throw an exception of type ExcUnitShapeValuesDoNotExist if the shape functions of the FiniteElement under consideration depend on the shape of the cell in real space, i.e., if the shape functions are not defined by mapping from the reference cell. Some non-conforming elements are defined this way, as is the FE_DGPNonparametric class, to name just one example.
The default implementation of this virtual function does exactly this, i.e., it simply throws an exception of type ExcUnitShapeValuesDoNotExist.
Reimplemented in FE_DGPNonparametric< dim, spacedim >, FE_NedelecSZ< dim, spacedim >, FE_Poly< dim, spacedim >, FE_Poly< dim, dim >, FE_Poly< dim, spacedim >, FE_PolyTensor< dim, spacedim >, FE_PolyTensor< dim, dim >, FE_PolyTensor< dim, spacedim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.
template<
int dim,
int spacedim = dim>
Return the tensor of third derivatives of the ith
shape function at point p
on the unit cell. The derivatives are derivatives on the unit cell with respect to unit cell coordinates. If the finite element is vector-valued, then return the value of the only non-zero component of the vector value of this shape function. If the shape function has more than one non-zero component (which we refer to with the term non-primitive), then derived classes implementing this function should throw an exception of type ExcShapeFunctionNotPrimitive. In that case, use the shape_3rd_derivative_component() function.
Implementations of this function should throw an exception of type ExcUnitShapeValuesDoNotExist if the shape functions of the FiniteElement under consideration depend on the shape of the cell in real space, i.e., if the shape functions are not defined by mapping from the reference cell. Some non-conforming elements are defined this way, as is the FE_DGPNonparametric class, to name just one example.
The default implementation of this virtual function does exactly this, i.e., it simply throws an exception of type ExcUnitShapeValuesDoNotExist.
Reimplemented in FE_Poly< dim, spacedim >, FE_Poly< dim, dim >, FE_Poly< dim, spacedim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.
template<
int dim,
int spacedim = dim>
Return the tensor of fourth derivatives of the ith
shape function at point p
on the unit cell. The derivatives are derivatives on the unit cell with respect to unit cell coordinates. If the finite element is vector-valued, then return the value of the only non-zero component of the vector value of this shape function. If the shape function has more than one non-zero component (which we refer to with the term non-primitive), then derived classes implementing this function should throw an exception of type ExcShapeFunctionNotPrimitive. In that case, use the shape_4th_derivative_component() function.
Implementations of this function should throw an exception of type ExcUnitShapeValuesDoNotExist if the shape functions of the FiniteElement under consideration depend on the shape of the cell in real space, i.e., if the shape functions are not defined by mapping from the reference cell. Some non-conforming elements are defined this way, as is the FE_DGPNonparametric class, to name just one example.
The default implementation of this virtual function does exactly this, i.e., it simply throws an exception of type ExcUnitShapeValuesDoNotExist.
Reimplemented in FE_Poly< dim, spacedim >, FE_Poly< dim, dim >, FE_Poly< dim, spacedim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.
template<
int dim,
int spacedim = dim>
Return the matrix that describes restricting a finite element field from the given child
(as obtained by the given refinement_case
) to the parent cell. The interpretation of the returned matrix depends on what restriction_is_additive() returns for each shape function.
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 an exception of type FiniteElement::ExcProjectionVoid. You can check whether this would happen by first calling the restriction_is_implemented() or the isotropic_restriction_is_implemented() function.
- Note
- The term "restriction" is also used in the definition of multigrid methods, but in a different context. For the current function, we are interested in the interpolation of a finite element function from child cells to the parent cell, and this results in a situation where you have a function that is constant on a fine mesh, the function "restricted" to the coarse mesh will also be constant and have the same value. In contrast, in the multigrid context, the restriction operation is applied to residual (rather than solution) vectors, and is usually chosen as the transpose of the prolongation operation. For the multigrid restriction operation, one would consequently not typically expect that it maps a constant function to a constant function. In other words, the meaning of the term "restriction" in the current context and in the multigrid context is quite different, and the two should not be confused.
Reimplemented in FE_Q_Bubbles< dim, spacedim >, FE_Bernstein< dim, spacedim >, FE_DGQ< dim, spacedim >, FE_DGQ< dim, dim >, FE_Enriched< dim, spacedim >, FE_Nedelec< dim >, FE_Q_Base< dim, spacedim >, FE_Q_Base< dim, dim >, FE_Q_Base< dim, spacedim >, FE_RaviartThomasNodal< dim >, FE_SimplexPoly< dim, spacedim >, FE_SimplexPoly< dim, dim >, FE_SimplexDGP< dim, spacedim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.
template<
int dim,
int spacedim = dim>
Prolongation/embedding matrix between grids.
The identity operator from a coarse grid space into a fine grid space (where both spaces are identified as functions defined on the parent and child cells) is associated with a matrix P
that maps the corresponding representations of these functions in terms of their nodal values. 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 prolongation matrices are not implemented in the derived finite element class, this function aborts with an exception of type FiniteElement::ExcEmbeddingVoid. You can check whether this would happen by first calling the prolongation_is_implemented() or the isotropic_prolongation_is_implemented() function.
Reimplemented in FE_Q_Bubbles< dim, spacedim >, FE_Bernstein< dim, spacedim >, FE_DGQ< dim, spacedim >, FE_DGQ< dim, dim >, FE_Enriched< dim, spacedim >, FE_Nedelec< dim >, FE_NedelecSZ< dim, spacedim >, FE_Q_Base< dim, spacedim >, FE_Q_Base< dim, dim >, FE_Q_Base< dim, spacedim >, FE_Q_Hierarchical< dim >, FE_RaviartThomasNodal< dim >, FE_SimplexPoly< dim, spacedim >, FE_SimplexPoly< dim, dim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.
template<
int dim,
int spacedim = dim>
Return the matrix interpolating from the given finite element to the present one. The size of the matrix is then dofs_per_cell times source.dofs_per_cell
.
Derived elements will have to implement this function. They may only provide interpolation matrices for certain source finite elements, for example those from the same family. If they don't implement interpolation from a given element, then they must throw an exception of type ExcInterpolationNotImplemented.
Reimplemented in FE_Bernstein< dim, spacedim >, FE_DGQ< dim, spacedim >, FE_DGQ< dim, dim >, FE_Q_Base< dim, spacedim >, FE_Q_Base< dim, dim >, FE_Q_Base< dim, spacedim >, FE_Q_Bubbles< dim, spacedim >, FE_Q_DG0< dim, spacedim >, FESystem< dim, spacedim >, FESystem< dim, dim >, FESystem< dim, spacedim >, FE_Nothing< dim, spacedim >, and FE_Nothing< dim, dim >.
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>
Given an index in the natural ordering of indices on a face, return the index of the same degree of freedom on the cell.
To explain the concept, consider the case where we would like to know whether a degree of freedom on a face, for example as part of an FESystem element, is primitive. Unfortunately, the is_primitive() function in the FiniteElement class takes a cell index, so we would need to find the cell index of the shape function that corresponds to the present face index. This function does that.
Code implementing this would then look like this:
if (fe.is_primitive(fe.face_to_cell_index(i,
some_face_no)))
const unsigned int dofs_per_face
The function takes additional arguments that account for the fact that actual faces can be in their standard ordering with respect to the cell under consideration, or can be flipped, oriented, etc.
- Parameters
-
face_dof_index | The index of the degree of freedom on a face. This index must be between zero and dofs_per_face. |
face | The number of the face this degree of freedom lives on. This number must be between zero and GeometryInfo::faces_per_cell. |
combined_orientation | The combined orientation flag containing the orientation, rotation, and flip of the face. See GlossCombinedOrientation. |
- Returns
- The index of this degree of freedom within the set of degrees of freedom on the entire cell. The returned value will be between zero and dofs_per_cell.
- Note
- This function exists in this class because that is where it was first implemented. However, it can't really work in the most general case without knowing what element we have. The reason is that when a face is flipped or rotated, we also need to know whether we need to swap the degrees of freedom on this face, or whether they are immune from this. For this, consider the situation of a \(Q_3\) element in 2d. If face_flip is true, then we need to consider the two degrees of freedom on the edge in reverse order. On the other hand, if the element were a \(Q_1^2\), then because the two degrees of freedom on this edge belong to different vector components, they should not be considered in reverse order. What all of this shows is that the function can't work if there are more than one degree of freedom per line or quad, and that in these cases the function will throw an exception pointing out that this functionality will need to be provided by a derived class that knows what degrees of freedom actually represent.
Reimplemented in FE_Q_Base< dim, spacedim >, FE_Q_Base< dim, dim >, FE_Q_Base< dim, spacedim >, FE_SimplexPoly< dim, spacedim >, FE_SimplexPoly< dim, dim >, FESystem< dim, spacedim >, FESystem< dim, dim >, and FESystem< dim, spacedim >.
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).