Classes | |
class | InternalData |
Public Member Functions | |
FE_PolyTensor (const unsigned int degree, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< std::vector< bool > > &nonzero_components) | |
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 unsigned int | n_base_elements () const |
virtual const FiniteElement < dim, spacedim > & | base_element (const unsigned int index) const |
virtual unsigned int | element_multiplicity (const unsigned int index) const |
virtual UpdateFlags | update_once (const UpdateFlags flags) const |
virtual UpdateFlags | update_each (const UpdateFlags flags) const |
Protected Member Functions | |
virtual Mapping< dim, spacedim > ::InternalDataBase * | get_data (const UpdateFlags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature) const |
virtual void | fill_fe_values (const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data, enum CellSimilarity::Similarity &cell_similarity) const |
virtual void | fill_fe_face_values (const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const |
virtual void | fill_fe_subface_values (const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const |
Protected Attributes | |
MappingType | mapping_type |
POLY | poly_space |
FullMatrix< double > | inverse_node_matrix |
Point< dim > | cached_point |
std::vector< Tensor< 1, dim > > | cached_values |
std::vector< Tensor< 2, dim > > | cached_grads |
std::vector< Tensor< 3, dim > > | cached_grad_grads |
This class gives a unified framework for the implementation of FiniteElement classes based on Tensor valued polynomial spaces like PolynomialsBDM and PolynomialsRaviartThomas.
Every class that implements following function can be used as template parameter POLY.
void compute (const Point<dim> &unit_point, std::vector<Tensor<1,dim> > &values, std::vector<Tensor<2,dim> > &grads, std::vector<Tensor<3,dim> > &grad_grads) const;
In many cases, the node functionals depend on the shape of the mesh cell, since they evaluate normal or tangential components on the faces. In order to allow for a set of transformations, the variable mapping_type has been introduced. It should also be set in the constructor of a derived class.
This class is not a fully implemented FiniteElement class, but implements some common features of vector valued elements based on vector valued polynomial classes. What's missing here in particular is information on the topological location of the node values.
For more information on the template parameter spacedim
, see the documentation for the class Triangulation.
Any derived class must decide on the polynomial space to use. This polynomial space should be implemented simply as a set of vector valued polynomials like PolynomialsBDM and PolynomialsRaviartThomas. In order to facilitate this implementation, the basis of this space may be arbitrary.
In most cases, the set of desired node values and the basis functions
will not fulfill the interpolation condition
.
The use of the member data inverse_node_matrix allows to compute the basis automatically, after the node values for each original basis function of the polynomial space have been computed.
Therefore, the constructor of a derived class should have a structure like this (example for interpolation in support points):
* fill_support_points(); * * const unsigned int n_dofs = this->dofs_per_cell; * FullMatrix<double> N(n_dofs, n_dofs); * * for (unsigned int i=0;i<n_dofs;++i) * { * const Point<dim>& p = this->unit_support_point[i]; * * for (unsigned int j=0;j<n_dofs;++j) * for (unsigned int d=0;d<dim;++d) * N(i,j) += node_vector[i][d] * * this->shape_value_component(j, p, d); * } * * this->inverse_node_matrix.reinit(n_dofs, n_dofs); * this->inverse_node_matrix.invert(N); *
In most cases, vector valued basis functions must be transformed when mapped from the reference cell to the actual grid cell. These transformations can be selected from the set MappingType and stored in mapping_type. Therefore, each constructor should contain a line like:
* this->mapping_type = this->mapping_none; *
FE_PolyTensor< POLY, dim, spacedim >::FE_PolyTensor | ( | const unsigned int | degree, | |
const FiniteElementData< dim > & | fe_data, | |||
const std::vector< bool > & | restriction_is_additive_flags, | |||
const std::vector< std::vector< bool > > & | nonzero_components | |||
) |
Constructor.
degree:
constructor argument for poly. May be different from fe_data.degree
. virtual double FE_PolyTensor< POLY, dim, spacedim >::shape_value | ( | const unsigned int | i, | |
const Point< dim > & | p | |||
) | const [virtual] |
Since these elements are vector valued, an exception is thrown.
Reimplemented from FiniteElement< dim, spacedim >.
virtual double FE_PolyTensor< POLY, dim, spacedim >::shape_value_component | ( | const unsigned int | i, | |
const Point< dim > & | p, | |||
const unsigned int | component | |||
) | const [virtual] |
Just like for shape_value(), but this function will be called when the shape function has more than one non-zero vector component. In that case, this function should return the value of the component-th
vector component of the ith
shape function at point p
.
Reimplemented from FiniteElement< dim, spacedim >.
virtual Tensor<1,dim> FE_PolyTensor< POLY, dim, spacedim >::shape_grad | ( | const unsigned int | i, | |
const Point< dim > & | p | |||
) | const [virtual] |
Since these elements are vector valued, an exception is thrown.
Reimplemented from FiniteElement< dim, spacedim >.
virtual Tensor<1,dim> FE_PolyTensor< POLY, dim, spacedim >::shape_grad_component | ( | const unsigned int | i, | |
const Point< dim > & | p, | |||
const unsigned int | component | |||
) | const [virtual] |
Just like for shape_grad(), but this function will be called when the shape function has more than one non-zero vector component. In that case, this function should return the gradient of the component-th
vector component of the ith
shape function at point p
.
Reimplemented from FiniteElement< dim, spacedim >.
virtual Tensor<2,dim> FE_PolyTensor< POLY, dim, spacedim >::shape_grad_grad | ( | const unsigned int | i, | |
const Point< dim > & | p | |||
) | const [virtual] |
Since these elements are vector valued, an exception is thrown.
Reimplemented from FiniteElement< dim, spacedim >.
virtual Tensor<2,dim> FE_PolyTensor< POLY, dim, spacedim >::shape_grad_grad_component | ( | const unsigned int | i, | |
const Point< dim > & | p, | |||
const unsigned int | component | |||
) | const [virtual] |
Just like for shape_grad_grad(), but this function will be called when the shape function has more than one non-zero vector component. In that case, this function should return the gradient of the component-th
vector component of the ith
shape function at point p
.
Reimplemented from FiniteElement< dim, spacedim >.
virtual unsigned int FE_PolyTensor< POLY, dim, spacedim >::n_base_elements | ( | ) | const [virtual] |
Number of base elements in a mixed discretization. Since this is not a composed element, return one.
Implements FiniteElement< dim, spacedim >.
Reimplemented in FE_ABF< dim >.
virtual const FiniteElement<dim,spacedim>& FE_PolyTensor< POLY, dim, spacedim >::base_element | ( | const unsigned int | index | ) | const [virtual] |
Access to base element objects. Since this element is not composed of several elements, base_element(0)
is this
, and all other indices throw an error.
Implements FiniteElement< dim, spacedim >.
Reimplemented in FE_ABF< dim >.
virtual unsigned int FE_PolyTensor< POLY, dim, spacedim >::element_multiplicity | ( | const unsigned int | index | ) | const [virtual] |
Multiplicity of base element index
. Since this is not a composed element, element_multiplicity(0)
returns one, and all other indices will throw an error.
Implements FiniteElement< dim, spacedim >.
Reimplemented in FE_ABF< dim >.
virtual UpdateFlags FE_PolyTensor< POLY, dim, spacedim >::update_once | ( | const UpdateFlags | flags | ) | const [virtual] |
Given flags
, determines the values which must be computed only for the reference cell. Make sure, that mapping_type is set by the derived class, such that this function can operate correctly.
Implements FiniteElement< dim, spacedim >.
Reimplemented in FE_ABF< dim >.
virtual UpdateFlags FE_PolyTensor< POLY, dim, spacedim >::update_each | ( | const UpdateFlags | flags | ) | const [virtual] |
Given flags
, determines the values which must be computed in each cell cell. Make sure, that mapping_type is set by the derived class, such that this function can operate correctly.
Implements FiniteElement< dim, spacedim >.
Reimplemented in FE_ABF< dim >.
virtual Mapping<dim,spacedim>::InternalDataBase* FE_PolyTensor< POLY, dim, spacedim >::get_data | ( | const | flags, | |
const Mapping< dim, spacedim > & | mapping, | |||
const Quadrature< dim > & | quadrature | |||
) | const [protected, virtual] |
Prepare internal data structures and fill in values independent of the cell. Returns a pointer to an object of which the caller of this function then has to assume ownership (which includes destruction when it is no more needed).
Implements FiniteElement< dim, spacedim >.
virtual void FE_PolyTensor< POLY, dim, spacedim >::fill_fe_values | ( | const Mapping< dim, spacedim > & | mapping, | |
const typename Triangulation< dim, spacedim >::cell_iterator & | cell, | |||
const Quadrature< dim > & | quadrature, | |||
typename Mapping< dim, spacedim >::InternalDataBase & | mapping_internal, | |||
typename Mapping< dim, spacedim >::InternalDataBase & | fe_internal, | |||
FEValuesData< dim, spacedim > & | data, | |||
enum CellSimilarity::Similarity & | cell_similarity | |||
) | const [protected, virtual] |
Fill the fields of FEValues. This function performs all the operations needed to compute the data of an FEValues object.
The same function in mapping
must have been called for the same cell first!
Implements FiniteElement< dim, spacedim >.
virtual void FE_PolyTensor< POLY, dim, spacedim >::fill_fe_face_values | ( | const Mapping< dim, spacedim > & | mapping, | |
const typename Triangulation< dim, spacedim >::cell_iterator & | cell, | |||
const unsigned int | face_no, | |||
const Quadrature< dim-1 > & | quadrature, | |||
typename Mapping< dim, spacedim >::InternalDataBase & | mapping_internal, | |||
typename Mapping< dim, spacedim >::InternalDataBase & | fe_internal, | |||
FEValuesData< dim, spacedim > & | data | |||
) | const [protected, virtual] |
Fill the fields of FEFaceValues. This function performs all the operations needed to compute the data of an FEFaceValues object.
The same function in mapping
must have been called for the same cell first!
Implements FiniteElement< dim, spacedim >.
virtual void FE_PolyTensor< POLY, dim, spacedim >::fill_fe_subface_values | ( | const Mapping< dim, spacedim > & | mapping, | |
const typename Triangulation< dim, spacedim >::cell_iterator & | cell, | |||
const unsigned int | face_no, | |||
const unsigned int | sub_no, | |||
const Quadrature< dim-1 > & | quadrature, | |||
typename Mapping< dim, spacedim >::InternalDataBase & | mapping_internal, | |||
typename Mapping< dim, spacedim >::InternalDataBase & | fe_internal, | |||
FEValuesData< dim, spacedim > & | data | |||
) | const [protected, virtual] |
Fill the fields of FESubfaceValues. This function performs all the operations needed to compute the data of an FESubfaceValues object.
The same function in mapping
must have been called for the same cell first!
Implements FiniteElement< dim, spacedim >.
MappingType FE_PolyTensor< POLY, dim, spacedim >::mapping_type [protected] |
The mapping type to be used to map shape functions from the reference cell to the mesh cell.
POLY FE_PolyTensor< POLY, dim, spacedim >::poly_space [protected] |
The polynomial space. Its type is given by the template parameter POLY.
FullMatrix<double> FE_PolyTensor< POLY, dim, spacedim >::inverse_node_matrix [protected] |
The inverse of the matrix aij of node values Ni applied to polynomial pj. This matrix is used to convert polynomials in the "raw" basis provided in poly_space to the basis dual to the node functionals on the reference cell.
This object is not filled by FE_PolyTensor, but is a chance for a derived class to allow for reorganization of the basis functions. If it is left empty, the basis in poly_space is used.
Point<dim> FE_PolyTensor< POLY, dim, spacedim >::cached_point [mutable, protected] |
If a shape function is computed at a single point, we must compute all of them to apply inverse_node_matrix. In order to avoid too much overhead, we cache the point and the function values for the next evaluation.
std::vector<Tensor<1,dim> > FE_PolyTensor< POLY, dim, spacedim >::cached_values [mutable, protected] |
Cached shape function values after call to shape_value_component().
std::vector<Tensor<2,dim> > FE_PolyTensor< POLY, dim, spacedim >::cached_grads [mutable, protected] |
Cached shape function gradients after call to shape_grad_component().
std::vector<Tensor<3,dim> > FE_PolyTensor< POLY, dim, spacedim >::cached_grad_grads [mutable, protected] |
Cached second derivatives of shape functions after call to shape_grad_grad_component().