![]() |
Reference documentation for deal.II version 8.1.0
|
#include <mapping_q.h>
Classes | |
class | InternalData |
Public Member Functions | |
MappingQ (const unsigned int p, const bool use_mapping_q_on_all_cells=false) | |
MappingQ (const MappingQ< dim, spacedim > &mapping) | |
virtual | ~MappingQ () |
virtual Point< spacedim > | transform_unit_to_real_cell (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const |
virtual Point< dim > | transform_real_to_unit_cell (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const |
virtual void | transform (const VectorSlice< const std::vector< Tensor< 1, dim > > > input, VectorSlice< std::vector< Tensor< 1, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const MappingType type) const |
virtual void | transform (const VectorSlice< const std::vector< DerivativeForm< 1, dim, spacedim > > > input, VectorSlice< std::vector< Tensor< 2, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const MappingType type) const |
virtual void | transform (const VectorSlice< const std::vector< Tensor< 2, dim > > > input, VectorSlice< std::vector< Tensor< 2, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const MappingType type) const |
unsigned int | get_degree () const |
virtual Mapping< dim, spacedim > * | clone () const |
![]() | |
MappingQ1 () | |
void | compute_shapes (const std::vector< Point< dim > > &unit_points, InternalData &data) const |
void | compute_data (const UpdateFlags flags, const Quadrature< dim > &quadrature, const unsigned int n_orig_q_points, InternalData &data) const |
void | compute_face_data (const UpdateFlags flags, const Quadrature< dim > &quadrature, const unsigned int n_orig_q_points, InternalData &data) const |
void | compute_fill (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int npts, const DataSetDescriptor data_set, const CellSimilarity::Similarity cell_similarity, InternalData &data, std::vector< Point< spacedim > > &quadrature_points) const |
void | compute_fill_face (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const unsigned int npts, const DataSetDescriptor data_set, const std::vector< double > &weights, InternalData &mapping_data, 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 |
virtual void | compute_shapes_virtual (const std::vector< Point< dim > > &unit_points, InternalData &data) const |
Point< spacedim > | transform_unit_to_real_cell_internal (const InternalData &mdata) const |
Point< dim > | transform_real_to_unit_cell_internal (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p, const Point< dim > &initial_p_unit, InternalData &mdata) const |
virtual bool | preserves_vertex_locations () const |
template<> | |
Point< 2 > | transform_real_to_unit_cell_internal (const Triangulation< 2, 3 >::cell_iterator &cell, const Point< 3 > &p, const Point< 2 > &initial_p_unit, InternalData &mdata) const |
template<> | |
Point< 1 > | transform_real_to_unit_cell_internal (const Triangulation< 1, 2 >::cell_iterator &cell, const Point< 2 > &p, const Point< 1 > &initial_p_unit, InternalData &mdata) const |
template<> | |
Point< 1 > | transform_real_to_unit_cell_internal (const Triangulation< 1, 3 >::cell_iterator &cell, const Point< 3 > &p, const Point< 1 > &initial_p_unit, InternalData &mdata) const |
![]() | |
virtual | ~Mapping () |
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 |
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) |
Protected Member Functions | |
virtual void | fill_fe_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename 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 |
virtual void | fill_fe_face_values (const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, typename std::vector< Tensor< 1, spacedim > > &exterior_form, typename std::vector< Point< spacedim > > &normal_vectors) 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, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, typename std::vector< Tensor< 1, spacedim > > &exterior_form, typename std::vector< Point< spacedim > > &normal_vectors) const |
virtual void | add_line_support_points (const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const |
virtual void | add_quad_support_points (const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const |
![]() | |
template<int rank> | |
void | transform_fields (const VectorSlice< const std::vector< Tensor< rank, dim > > > input, VectorSlice< std::vector< Tensor< rank, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const MappingType type) const |
template<int rank> | |
void | transform_gradients (const VectorSlice< const std::vector< Tensor< rank, dim > > > input, VectorSlice< std::vector< Tensor< rank, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const MappingType type) const |
template<int rank> | |
void | transform_differential_forms (const VectorSlice< const std::vector< DerivativeForm< rank, dim, spacedim > > > input, VectorSlice< std::vector< DerivativeForm< rank, spacedim, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const MappingType type) const |
template<int dim_> | |
Point< dim_ > | transform_real_to_unit_cell_internal_codim1 (const typename Triangulation< dim_, dim_+1 >::cell_iterator &cell, const Point< dim_+1 > &p, const Point< dim_ > &initial_p_unit, InternalData &mdata) const |
Point< dim > | transform_real_to_unit_cell_initial_guess (const std::vector< Point< spacedim > > &vertex, const Point< spacedim > &p) const |
Private Member Functions | |
virtual Mapping< dim, spacedim >::InternalDataBase * | get_data (const UpdateFlags, const Quadrature< dim > &quadrature) const |
virtual Mapping< dim, spacedim >::InternalDataBase * | get_face_data (const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const |
virtual Mapping< dim, spacedim >::InternalDataBase * | get_subface_data (const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const |
virtual void | compute_shapes_virtual (const std::vector< Point< dim > > &unit_points, typename MappingQ1< dim, spacedim >::InternalData &data) const |
void | set_laplace_on_quad_vector (Table< 2, double > &loqvs) const |
void | set_laplace_on_hex_vector (Table< 2, double > &lohvs) const |
void | compute_laplace_vector (Table< 2, double > &lvs) const |
void | apply_laplace_vector (const Table< 2, double > &lvs, std::vector< Point< spacedim > > &a) const |
virtual void | compute_mapping_support_points (const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const |
void | compute_support_points_laplace (const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const |
DeclException1 (ExcLaplaceVectorNotSet, int,<< "laplace_vector not set for degree="<< arg1<< ".") | |
Private Attributes | |
Table< 2, double > | laplace_on_quad_vector |
Table< 2, double > | laplace_on_hex_vector |
const unsigned int | degree |
const unsigned int | n_inner |
const unsigned int | n_outer |
const TensorProductPolynomials< dim > * | tensor_pols |
const unsigned int | n_shape_functions |
const std::vector< unsigned int > | renumber |
const bool | use_mapping_q_on_all_cells |
const FE_Q< dim > | feq |
Friends | |
template<int , int > | |
class | MappingQ |
Additional Inherited Members | |
![]() | |
typedef QProjector< dim >::DataSetDescriptor | DataSetDescriptor |
Mapping class that uses Qp-mappings on boundary cells. The mapping shape functions make use of tensor product polynomials with unit cell support points equal to the points of the Gauss-Lobatto quadrature formula. These points give a well-conditioned interpolation also for very high orders and are therefore preferred over equidistant support points.
For more details about Qp-mappings, see the `mapping' report at deal.II/doc/reports/mapping_q/index.html
in the `Reports' section of `Documentation'.
For more information about the spacedim
template parameter check the documentation of FiniteElement or the one of Triangulation.
Definition at line 28 of file fe_dgp_nonparametric.h.
MappingQ< dim, spacedim >::MappingQ | ( | const unsigned int | p, |
const bool | use_mapping_q_on_all_cells = false |
||
) |
Constructor. p
gives the degree of mapping polynomials on boundary cells.
The second argument determines whether the higher order mapping should also be used on interior cells. If its value is false
(the default), the a lower-order mapping is used in the interior. This is sufficient for most cases where higher order mappings are only used to better approximate the boundary. In that case, cells bounded by straight lines are acceptable in the interior. However, there are cases where one would also like to use a higher order mapping in the interior. The MappingQEulerian class is one such case.
MappingQ< dim, spacedim >::MappingQ | ( | const MappingQ< dim, spacedim > & | mapping | ) |
Copy constructor. Performs a deep copy, i.e. duplicates what tensor_pols points to instead of simply copying the tensor_pols pointer as done by a default copy constructor.
Destructor.
|
virtual |
Transforms the point p
on the unit cell to the point p_real
on the real cell cell
and returns p_real
.
Reimplemented from MappingQ1< dim, spacedim >.
|
virtual |
Transforms the point p
on the real cell to the point p_unit
on the unit cell cell
and returns p_unit
.
Uses Newton iteration and the transform_unit_to_real_cell
function.
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. Reimplemented from MappingQ1< dim, spacedim >.
Return the degree of the mapping, i.e. the value which was passed to the constructor.
|
virtual |
Return a pointer to a copy of the present object. The caller of this copy then assumes ownership of it.
Reimplemented from MappingQ1< dim, spacedim >.
Reimplemented in MappingQEulerian< dim, VECTOR, spacedim >, and MappingC1< dim, spacedim >.
|
protectedvirtual |
Implementation of the interface in Mapping.
Reimplemented from MappingQ1< dim, spacedim >.
Reimplemented in MappingQEulerian< dim, VECTOR, spacedim >.
|
protectedvirtual |
Implementation of the interface in Mapping.
Reimplemented from MappingQ1< dim, spacedim >.
|
protectedvirtual |
Implementation of the interface in Mapping.
Reimplemented from MappingQ1< dim, spacedim >.
|
protectedvirtual |
For dim=2,3
. Append the support points of all shape functions located on bounding lines to the vector a
. Points located on the line but not on vertices are not included.
Needed by the compute_support_points_laplace
function . For dim=1
this function is empty.
This function is made virtual in order to allow derived classes to choose shape function support points differently than the present class, which chooses the points as interpolation points on the boundary.
|
protectedvirtual |
For dim=3
. Append the support points of all shape functions located on bounding faces (quads in 3d) to the vector a
. Points located on the quad but not on vertices are not included.
Needed by the compute_support_points_laplace
function. For dim=1
and dim=2
this function is empty.
This function is made virtual in order to allow derived classes to choose shape function support points differently than the present class, which chooses the points as interpolation points on the boundary.
|
privatevirtual |
Prepare internal data structures and fill in values independent of the cell.
Reimplemented from MappingQ1< dim, spacedim >.
|
privatevirtual |
Prepare internal data structure for transformation of faces and fill in values independent of the cell.
Reimplemented from MappingQ1< dim, spacedim >.
|
privatevirtual |
Prepare internal data structure for transformation of children of faces and fill in values independent of the cell.
Reimplemented from MappingQ1< dim, spacedim >.
|
privatevirtual |
Compute shape values and/or derivatives.
|
private |
This function is needed by the constructor of MappingQ<dim,spacedim>
for dim=
2 and 3.
For degree<4
this function sets the laplace_on_quad_vector
to the hardcoded data. For degree>=4
and MappingQ<2> this vector is computed.
For the definition of the laplace_on_quad_vector
please refer to equation (8) of the `mapping' report.
|
private |
This function is needed by the constructor of MappingQ<3>
.
For degree==2
this function sets the laplace_on_hex_vector
to the hardcoded data. For degree>2
this vector is computed.
For the definition of the laplace_on_hex_vector
please refer to equation (8) of the `mapping' report.
|
private |
Computes the laplace_on_quad(hex)_vector
.
Called by the set_laplace_on_quad(hex)_vector
functions if the data is not yet hardcoded.
For the definition of the laplace_on_quad(hex)_vector
please refer to equation (8) of the `mapping' report.
|
private |
Takes a laplace_on_hex(quad)_vector
and applies it to the vector a
to compute the inner support points as a linear combination of the exterior points.
The vector a
initially contains the locations of the n_outer
points, the n_inner
computed inner points are appended.
See equation (7) of the `mapping' report.
|
privatevirtual |
Computes the support points of the mapping.
Reimplemented from MappingQ1< dim, spacedim >.
Reimplemented in MappingQEulerian< dim, VECTOR, spacedim >.
|
private |
Computes all support points of the mapping shape functions. The inner support points (ie. support points in quads for 2d, in hexes for 3d) are computed using the solution of a Laplace equation with the position of the outer support points as boundary values, in order to make the transformation as smooth as possible.
|
private |
Exception.
Declare other MappingQ classes friends.
Definition at line 440 of file mapping_q.h.
|
private |
Needed by the laplace_on_quad
function (for dim==2
). Filled by the constructor.
Sizes: laplace_on_quad_vector.size()= number of inner unit_support_points laplace_on_quad_vector[i].size()= number of outer unit_support_points, i.e. unit_support_points on the boundary of the quad
For the definition of this vector see equation (8) of the `mapping' report.
Definition at line 371 of file mapping_q.h.
|
private |
Needed by the laplace_on_hex
function (for dim==3
). Filled by the constructor.
For the definition of this vector see equation (8) of the `mapping' report.
Definition at line 380 of file mapping_q.h.
Degree p
of the polynomials used as shape functions for the Qp mapping of cells at the boundary.
Definition at line 393 of file mapping_q.h.
Number of inner mapping shape functions.
Definition at line 398 of file mapping_q.h.
Number of mapping shape functions on the boundary.
Definition at line 403 of file mapping_q.h.
|
private |
Pointer to the dim-dimensional
tensor product polynomials used as shape functions for the Qp mapping of cells at the boundary.
Definition at line 409 of file mapping_q.h.
|
private |
Number of the Qp tensor product shape functions.
Definition at line 414 of file mapping_q.h.
|
private |
Mapping from lexicographic to to the Qp shape function numbering. Its size is dofs_per_cell
.
Definition at line 420 of file mapping_q.h.
|
private |
If this flag is set true
then MappingQ
is used on all cells, not only on boundary cells.
Definition at line 426 of file mapping_q.h.
An FE_Q object which is only needed in 3D, since it knows how to reorder shape functions/DoFs on non-standard faces. This is used to reorder support points in the same way. We could make this a pointer to prevent construction in 1D and 2D, but since memory and time requirements are not particularly high this seems unnecessary at the moment.
Definition at line 435 of file mapping_q.h.