
Classes | |
| class | InternalDataBase |
Public Member Functions | |
| virtual | ~Mapping () |
| virtual Point< spacedim > | transform_unit_to_real_cell (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const =0 |
| virtual Point< dim > | transform_real_to_unit_cell (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const =0 |
| virtual void | transform (const VectorSlice< const std::vector< Tensor< 1, dim > > > input, VectorSlice< std::vector< Tensor< 1, spacedim > > > output, const InternalDataBase &internal, const MappingType type) const =0 |
| virtual void | transform (const VectorSlice< const std::vector< Tensor< 2, dim > > > input, VectorSlice< std::vector< Tensor< 2, spacedim > > > output, const InternalDataBase &internal, const MappingType type) const =0 |
| void | transform_covariant (const VectorSlice< const std::vector< Tensor< 1, dim > > > input, const unsigned int offset, VectorSlice< std::vector< Tensor< 1, spacedim > > > output, const InternalDataBase &internal) const |
| void | transform_covariant (const VectorSlice< const std::vector< Tensor< 2, dim > > > input, const unsigned int offset, VectorSlice< std::vector< Tensor< 2, spacedim > > > output, const InternalDataBase &internal) const |
| void | transform_contravariant (const VectorSlice< const std::vector< Tensor< 1, dim > > > input, const unsigned int offset, VectorSlice< std::vector< Tensor< 1, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const |
| void | transform_contravariant (const VectorSlice< const std::vector< Tensor< 2, dim > > > intput, const unsigned int offset, const VectorSlice< std::vector< Tensor< 2, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const |
| const Point< spacedim > & | support_point_value (const unsigned int index, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const |
| const Tensor< 2, spacedim > & | support_point_gradient (const unsigned int index, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const |
| const Tensor< 2, spacedim > & | support_point_inverse_gradient (const unsigned int index, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const |
| virtual Mapping< dim, spacedim > * | clone () const =0 |
| virtual bool | preserves_vertex_locations () const =0 |
Static Public Member Functions | |
| ::ExceptionBase & | ExcInvalidData () |
Private Member Functions | |
| virtual UpdateFlags | update_once (const UpdateFlags) const =0 |
| virtual UpdateFlags | update_each (const UpdateFlags) const =0 |
| virtual InternalDataBase * | get_data (const UpdateFlags, const Quadrature< dim > &quadrature) const =0 |
| virtual InternalDataBase * | get_face_data (const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const =0 |
| virtual InternalDataBase * | get_subface_data (const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const =0 |
| virtual void | fill_fe_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, InternalDataBase &internal, std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, std::vector< Tensor< 2, spacedim > > &jacobians, std::vector< Tensor< 3, spacedim > > &jacobian_grads, std::vector< Tensor< 2, spacedim > > &inverse_jacobians, std::vector< Point< spacedim > > &cell_normal_vectors, CellSimilarity::Similarity &cell_similarity) const =0 |
| virtual void | fill_fe_face_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, InternalDataBase &internal, std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, std::vector< Tensor< 1, spacedim > > &boundary_form, std::vector< Point< spacedim > > &normal_vectors) const =0 |
| 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, InternalDataBase &internal, std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, std::vector< Tensor< 1, spacedim > > &boundary_form, std::vector< Point< spacedim > > &normal_vectors) const =0 |
Friends | |
| class | FEValuesBase< dim, spacedim > |
| class | FEValues< dim, spacedim > |
| class | FEFaceValues< dim, spacedim > |
| class | FESubfaceValues< dim, spacedim > |
Abstract base class for mapping classes.
The interface for filling the tables of FEValues is provided. Everything else has to happen in derived classes.
The following paragraph applies to the implementation of FEValues. Usage of the class is as follows: first, call the functions update_once and update_each with the update flags you need. This includes the flags needed by the FiniteElement. Then call get_*_data and with the or'd results. This will initialize and return some internal data structures. On the first cell, call fill_fe_*_values with the result of update_once. Finally, on each cell, use fill_fe_*_values with the result of update_each to compute values for a special cell.
The mapping is a transformation
which maps the reference cell [0,1]dim to the actual grid cell. In order to describe the application of the mapping to different objects, we introduce the notation for the Jacobian
. For instance, in two dimensions, we have
Functions are simply mapped such that
Since finite element shape functions are usually defined on the reference cell, nothing needs to be done for them. For a function defined on the computational domain, the quadrature points need to be mapped, which is done in fill_fe_values() if update_quadrature_points is set in the update flags. The mapped quadrature points are then accessed through FEValuesBase::quadrature_point().
transform_quadrature_points for this.The volume form
is mapped such that for a grid cell Z
The transformed quadrature weights
are accessed through FEValuesBase::JxW() and computed in fill_fe_values(), if update_JxW_values is set in the update flags.
transform_quadrature_weights for this.Mappings of a vector field v and a tensor field T follows the general form
respectively, where the tensors A and B are determined by the MappingType enumerator. These transformations are performed through the two functions transform(). See the documentation there for possible choices.
A hint to implementators: no function except the two functions update_once and update_each may add any flags.
For more information about the spacedim template parameter check the documentation of FiniteElement or the one of Triangulation.
A general publication on differential geometry and finite elements is the survey
The description of the Piola transform has been taken from the lecture notes by Ronald H. W. Hoppe, University of Houston, Chapter 7.
Definition at line 223 of file mapping.h.
Virtual destructor.
| virtual Point<spacedim> Mapping< dim, spacedim >::transform_unit_to_real_cell | ( | const typename Triangulation< dim, spacedim >::cell_iterator & | cell, |
| const Point< dim > & | p | ||
| ) | const [pure virtual] |
Transforms the point p on the unit cell to the point p_real on the real cell cell and returns p_real.
Implemented in MappingCartesian< dim, spacedim >, MappingQ< dim, spacedim >, and MappingQ1< dim, spacedim >.
| virtual Point<dim> Mapping< dim, spacedim >::transform_real_to_unit_cell | ( | const typename Triangulation< dim, spacedim >::cell_iterator & | cell, |
| const Point< spacedim > & | p | ||
| ) | const [pure virtual] |
Transforms the point p on the real cell to the point p_unit on the unit cell cell and returns p_unit.
In the codimension one case, this function returns the normal projection of the real point p on the curve or surface identified by the cell.
Implemented in MappingCartesian< dim, spacedim >, MappingQ< dim, spacedim >, and MappingQ1< dim, spacedim >.
| virtual void Mapping< dim, spacedim >::transform | ( | const VectorSlice< const std::vector< Tensor< 1, dim > > > | input, |
| VectorSlice< std::vector< Tensor< 1, spacedim > > > | output, | ||
| const InternalDataBase & | internal, | ||
| const MappingType | type | ||
| ) | const [pure virtual] |
Transform a field of vectors or forms according to the selected MappingType.
The mapping types currently implemented by derived classes are:
A vector field on the reference cell is mapped to the physical cell by the Jacobian:
In physics, this is usually refered to as the contravariant transformation, and here it is triggered by the MappingType mapping_contravariant. Mathematically, it is the push forward of a vector field.
A field of one-forms on the reference cell is usually represented as a vector field as well, but is transformed by the pull back, which is selected by mapping_covariant:
Gradients of differentiable functions are transformed this way.
A field of n-1-forms on the reference cell is also represented by a vector field, but again transforms differently, namely by the Piola transform (mapping_piola)
| virtual void Mapping< dim, spacedim >::transform | ( | const VectorSlice< const std::vector< Tensor< 2, dim > > > | input, |
| VectorSlice< std::vector< Tensor< 2, spacedim > > > | output, | ||
| const InternalDataBase & | internal, | ||
| const MappingType | type | ||
| ) | const [pure virtual] |
Transform a field of rank two tensors according to the selected MappingType.
This function is most of the time applied to gradients of a Tensor<1,dim> object, thus in the formulas below, it is useful to think of
and
.
The mapping types currently implemented by derived classes are:
mapping_contravariant_gradient, which is the consistent gradient for a vector field mapped with mapping_contravariant, namely
Correspondingly, mapping_covariant_gradient is the consistent gradient of a one-form, namely,
| void Mapping< dim, spacedim >::transform_covariant | ( | const VectorSlice< const std::vector< Tensor< 1, dim > > > | input, |
| const unsigned int | offset, | ||
| VectorSlice< std::vector< Tensor< 1, spacedim > > > | output, | ||
| const InternalDataBase & | internal | ||
| ) | const |
| void Mapping< dim, spacedim >::transform_covariant | ( | const VectorSlice< const std::vector< Tensor< 2, dim > > > | input, |
| const unsigned int | offset, | ||
| VectorSlice< std::vector< Tensor< 2, spacedim > > > | output, | ||
| const InternalDataBase & | internal | ||
| ) | const |
| void Mapping< dim, spacedim >::transform_contravariant | ( | const VectorSlice< const std::vector< Tensor< 1, dim > > > | input, |
| const unsigned int | offset, | ||
| VectorSlice< std::vector< Tensor< 1, spacedim > > > | output, | ||
| const typename Mapping< dim, spacedim >::InternalDataBase & | internal | ||
| ) | const |
| void Mapping< dim, spacedim >::transform_contravariant | ( | const VectorSlice< const std::vector< Tensor< 2, dim > > > | intput, |
| const unsigned int | offset, | ||
| const VectorSlice< std::vector< Tensor< 2, spacedim > > > | output, | ||
| const typename Mapping< dim, spacedim >::InternalDataBase & | internal | ||
| ) | const |
| const Point<spacedim>& Mapping< dim, spacedim >::support_point_value | ( | const unsigned int | index, |
| const typename Mapping< dim, spacedim >::InternalDataBase & | internal | ||
| ) | const |
The transformed (generalized) support point.
| const Tensor<2,spacedim>& Mapping< dim, spacedim >::support_point_gradient | ( | const unsigned int | index, |
| const typename Mapping< dim, spacedim >::InternalDataBase & | internal | ||
| ) | const |
The Jacobian matrix of the transformation in the (generalized) support point.
| const Tensor<2,spacedim>& Mapping< dim, spacedim >::support_point_inverse_gradient | ( | const unsigned int | index, |
| const typename Mapping< dim, spacedim >::InternalDataBase & | internal | ||
| ) | const |
The inverse Jacobian matrix of the transformation in the (generalized) support point.
| virtual Mapping<dim,spacedim>* Mapping< dim, spacedim >::clone | ( | ) | const [pure virtual] |
Return a pointer to a copy of the present object. The caller of this copy then assumes ownership of it.
Since one can't create objects of class Mapping, this function of course has to be implemented by derived classes.
This function is mainly used by the hp::MappingCollection class.
Implemented in MappingC1< dim, spacedim >, MappingCartesian< dim, spacedim >, MappingQ< dim, spacedim >, MappingQ1< dim, spacedim >, MappingQ1Eulerian< dim, VECTOR, spacedim >, and MappingQEulerian< dim, VECTOR, spacedim >.
| virtual bool Mapping< dim, spacedim >::preserves_vertex_locations | ( | ) | const [pure virtual] |
Returns whether the mapping preserves vertex locations, i.e. whether the mapped location of the reference cell vertices (given by GeometryInfo::unit_cell_vertex()) equals the result of cell->vertex().
For example, implementations in derived classes return true for MappingQ, MappingQ1, MappingCartesian, but false for MappingQEulerian, MappingQ1Eulerian.
Implemented in MappingCartesian< dim, spacedim >, MappingQ1< dim, spacedim >, MappingQ1Eulerian< dim, VECTOR, spacedim >, and MappingQEulerian< dim, VECTOR, spacedim >.
| ::ExceptionBase& Mapping< dim, spacedim >::ExcInvalidData | ( | ) | [static] |
Exception
| virtual UpdateFlags Mapping< dim, spacedim >::update_once | ( | const | UpdateFlags ) | const [private, pure virtual] |
Indicate fields to be updated in the constructor of FEValues. Especially, fields not asked for by FEValues, but computed for efficiency reasons will be notified here.
See The interplay of UpdateFlags, Mapping and FiniteElement in FEValues.
Implemented in MappingCartesian< dim, spacedim >, and MappingQ1< dim, spacedim >.
| virtual UpdateFlags Mapping< dim, spacedim >::update_each | ( | const | UpdateFlags ) | const [private, pure virtual] |
The same as update_once(), but for the flags to be updated for each grid cell.
See The interplay of UpdateFlags, Mapping and FiniteElement in FEValues.
Implemented in MappingCartesian< dim, spacedim >, and MappingQ1< dim, spacedim >.
| virtual InternalDataBase* Mapping< dim, spacedim >::get_data | ( | const | UpdateFlags, |
| const Quadrature< dim > & | quadrature | ||
| ) | const [private, pure virtual] |
Prepare internal data structures and fill in values independent of the cell.
Implemented in MappingCartesian< dim, spacedim >, MappingQ< dim, spacedim >, and MappingQ1< dim, spacedim >.
| virtual InternalDataBase* Mapping< dim, spacedim >::get_face_data | ( | const UpdateFlags | flags, |
| const Quadrature< dim-1 > & | quadrature | ||
| ) | const [private, pure virtual] |
Prepare internal data structure for transformation of faces and fill in values independent of the cell.
Implemented in MappingCartesian< dim, spacedim >, MappingQ< dim, spacedim >, and MappingQ1< dim, spacedim >.
| virtual InternalDataBase* Mapping< dim, spacedim >::get_subface_data | ( | const UpdateFlags | flags, |
| const Quadrature< dim-1 > & | quadrature | ||
| ) | const [private, pure virtual] |
Prepare internal data structure for transformation of children of faces and fill in values independent of the cell.
Implemented in MappingCartesian< dim, spacedim >, MappingQ< dim, spacedim >, and MappingQ1< dim, spacedim >.
| virtual void Mapping< dim, spacedim >::fill_fe_values | ( | const typename Triangulation< dim, spacedim >::cell_iterator & | cell, |
| const Quadrature< dim > & | quadrature, | ||
| InternalDataBase & | internal, | ||
| std::vector< Point< spacedim > > & | quadrature_points, | ||
| std::vector< double > & | JxW_values, | ||
| std::vector< Tensor< 2, spacedim > > & | jacobians, | ||
| std::vector< Tensor< 3, spacedim > > & | jacobian_grads, | ||
| std::vector< Tensor< 2, spacedim > > & | inverse_jacobians, | ||
| std::vector< Point< spacedim > > & | cell_normal_vectors, | ||
| CellSimilarity::Similarity & | cell_similarity | ||
| ) | const [private, pure virtual] |
Fill the transformation fields of FEValues. Given a grid cell and the quadrature points on the unit cell, it computes all values specified by flags. The arrays to be filled have to have the correct size.
Values are split into two groups: first, quadrature_points and JxW_values are filled with the quadrature rule transformed to the cell in physical space.
The second group contains the matrices needed to transform vector-valued functions, namely jacobians, the derivatives jacobian_grads, and the inverse operations in inverse_jacobians. The function above adjusted with the variable cell_normal_vectors for the case of codimension 1
| virtual void Mapping< dim, spacedim >::fill_fe_face_values | ( | const typename Triangulation< dim, spacedim >::cell_iterator & | cell, |
| const unsigned int | face_no, | ||
| const Quadrature< dim-1 > & | quadrature, | ||
| InternalDataBase & | internal, | ||
| std::vector< Point< spacedim > > & | quadrature_points, | ||
| std::vector< double > & | JxW_values, | ||
| std::vector< Tensor< 1, spacedim > > & | boundary_form, | ||
| std::vector< Point< spacedim > > & | normal_vectors | ||
| ) | const [private, pure virtual] |
Performs the same as fill_fe_values on a face. Additionally, boundary_form (see GlossBoundaryForm) and normal_vectors can be computed on surfaces. Since the boundary form already contains the determinant of the Jacobian of the transformation, it is sometimes more economic to use the boundary form instead of the product of the unit normal and the transformed quadrature weight.
| virtual void Mapping< dim, spacedim >::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, | ||
| InternalDataBase & | internal, | ||
| std::vector< Point< spacedim > > & | quadrature_points, | ||
| std::vector< double > & | JxW_values, | ||
| std::vector< Tensor< 1, spacedim > > & | boundary_form, | ||
| std::vector< Point< spacedim > > & | normal_vectors | ||
| ) | const [private, pure virtual] |
See above.
friend class FEValuesBase< dim, spacedim > [friend] |
friend class FEValues< dim, spacedim > [friend] |
friend class FEFaceValues< dim, spacedim > [friend] |
friend class FESubfaceValues< dim, spacedim > [friend] |
documentation generated on Fri Feb 3 2012 06:04:09 by
doxygen
1.7.2