![]() |
Reference documentation for deal.II version 8.1.0
|
#include <mapping.h>
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< DerivativeForm< 1, dim, spacedim > > > input, VectorSlice< std::vector< Tensor< 2, 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 DEAL_II_DEPRECATED |
void | transform_covariant (const VectorSlice< const std::vector< DerivativeForm< 1, dim, spacedim > > > input, const unsigned int offset, VectorSlice< std::vector< Tensor< 2, spacedim > > > output, const InternalDataBase &internal) const DEAL_II_DEPRECATED |
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 DEAL_II_DEPRECATED |
void | transform_contravariant (const VectorSlice< const std::vector< DerivativeForm< 1, dim, spacedim > > > input, const unsigned int offset, const VectorSlice< std::vector< Tensor< 2, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const DEAL_II_DEPRECATED |
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 |
DeclException0 (ExcInvalidData) | |
DeclException0 (ExcTransformationFailed) | |
DeclException3 (ExcDistortedMappedCell, Point< spacedim >, double, int,<< "The image of the mapping applied to cell with center ["<< arg1<< "] is distorted. The cell geometry or the "<< "mapping are invalid, giving a non-positive volume "<< "fraction of "<< arg2<< " in quadrature point "<< arg3<< ".") | |
![]() | |
Subscriptor () | |
Subscriptor (const Subscriptor &) | |
virtual | ~Subscriptor () |
Subscriptor & | operator= (const Subscriptor &) |
void | subscribe (const char *identifier=0) const |
void | unsubscribe (const char *identifier=0) const |
unsigned int | n_subscriptions () const |
void | list_subscribers () const |
DeclException3 (ExcInUse, int, char *, std::string &,<< "Object of class "<< arg2<< " is still used by "<< arg1<< " other objects.\n"<< "(Additional information: "<< arg3<< ")\n"<< "Note the entry in the Frequently Asked Questions of "<< "deal.II (linked to from http://www.dealii.org/) for "<< "more information on what this error means.") | |
DeclException2 (ExcNoSubscriber, char *, char *,<< "No subscriber with identifier \""<< arg2<< "\" did subscribe to this object of class "<< arg1) | |
template<class Archive > | |
void | serialize (Archive &ar, const unsigned int version) |
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< DerivativeForm< 1, dim, spacedim > > &jacobians, std::vector< DerivativeForm< 2, dim, spacedim > > &jacobian_grads, std::vector< DerivativeForm< 1, spacedim, dim > > &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 Rspacedim. In order to describe the application of the mapping to different objects, we introduce the notation for the Jacobian
. For instance, if dim=spacedim=2, 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.The transfomation of vector fields, differential forms (gradients/jacobians) and gradients of vector fields between the reference cell and the actual grid cell follows the general form
where v is a vector field or a differential form and and T a tensor field of gradients. The differential forms A and B are determined by the MappingType enumerator. These transformations are performed through the 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 34 of file grid_out.h.
Virtual destructor.
|
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 >.
|
pure virtual |
Transforms the point p
on the real cell
to the corresponding point on the unit cell, and return its coordinates.
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
.
p
. If this is the case then this function throws an exception of type Mapping::ExcTransformationFailed . Whether the given point p
lies outside the cell can therefore be determined by checking whether the return reference coordinates lie inside of outside the reference cell (e.g., using GeometryInfo::is_inside_unit_cell) or whether the exception mentioned above has been thrown. Implemented in MappingCartesian< dim, spacedim >, MappingQ< dim, spacedim >, and MappingQ1< dim, spacedim >.
|
pure virtual |
Transform a field of vectors or 1-differential forms according to the selected MappingType.
mapping_bdm
, mapping_nedelec
, etc. This alias should be preferred to using the types below.The mapping types currently implemented by derived classes are:
mapping_contravariant:
maps a vector field on the reference cell is to the physical cell through the Jacobian:
In physics, this is usually referred to as the contravariant transformation. Mathematically, it is the push forward of a vector field.
mapping_covariant:
maps a field of one-forms on the reference cell to a field of one-forms on the physical cell. (theoretically this would refer to a DerivativeForm<1, dim, 1> but it canonically identified with a Tensor<1,dim>). Mathematically, it is the pull back of the differential form
In the case when dim=spacedim the previous formula reduces to
Gradients of scalar differentiable functions are transformed this way.
mapping_piola:
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
|
pure virtual |
Transform a field of differential forms from the reference cell to the physical cell.
It is useful to think of and
, with
a vector field.
The mapping types currently implemented by derived classes are:
mapping_covariant:
maps a field of forms on the reference cell to a field of forms on the physical cell. Mathematically, it is the pull back of the differential form
DerivativeForm<1, dim, rank>
. Unfortunately C++ does not allow templatized virtual functions. This is why we identify DerivativeForm<1, dim, 1>
with a Tensor<1,dim>
when using mapping_covariant() in the function transform above this one.
|
pure virtual |
Transform a tensor field from the reference cell to the physical cell. This tensors are most of times the jacobians in the reference cell of vector fields that have been pulled back from the physical cell.
The mapping types currently implemented by derived classes are:
mapping_contravariant_gradient
, it assumes so that
mapping_covariant_gradient
, it assumes so that
mapping_piola_gradient
, it assumes
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< DerivativeForm< 1, dim, spacedim > > > | 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< DerivativeForm< 1, dim, spacedim > > > | input, |
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.
|
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 MappingQ1< dim, spacedim >, MappingQ< dim, spacedim >, MappingCartesian< dim, spacedim >, MappingQ1Eulerian< dim, VECTOR, spacedim >, MappingQEulerian< dim, VECTOR, spacedim >, and MappingC1< dim, spacedim >.
|
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 MappingQ1< dim, spacedim >, MappingCartesian< dim, spacedim >, MappingQ1Eulerian< dim, VECTOR, spacedim >, and MappingQEulerian< dim, VECTOR, spacedim >.
Mapping< dim, spacedim >::DeclException0 | ( | ExcInvalidData | ) |
Exception
|
privatepure 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 MappingQ1< dim, spacedim >, and MappingCartesian< dim, spacedim >.
|
privatepure 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 MappingQ1< dim, spacedim >, and MappingCartesian< dim, spacedim >.
|
privatepure virtual |
Prepare internal data structures and fill in values independent of the cell.
Implemented in MappingQ1< dim, spacedim >, MappingQ< dim, spacedim >, and MappingCartesian< dim, spacedim >.
|
privatepure virtual |
Prepare internal data structure for transformation of faces and fill in values independent of the cell.
Implemented in MappingQ1< dim, spacedim >, MappingQ< dim, spacedim >, and MappingCartesian< dim, spacedim >.
|
privatepure virtual |
Prepare internal data structure for transformation of children of faces and fill in values independent of the cell.
Implemented in MappingQ1< dim, spacedim >, MappingQ< dim, spacedim >, and MappingCartesian< dim, spacedim >.
|
privatepure 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
|
privatepure 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.
|
privatepure virtual |
See above.
|
friend |