![]() |
Reference documentation for deal.II version GIT 3c37cfefb2 2023-03-24 04:00:03+00:00
|
#include <deal.II/fe/mapping_q.h>
Public Types | |
using | VectorizedArrayType = VectorizedArray< double, std::min< std::size_t >(VectorizedArray< double >::size(),(dim<=2 ? 2 :4))> |
Public Member Functions | |
InternalData (const unsigned int polynomial_degree) | |
void | initialize (const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points) |
void | initialize_face (const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points) |
void | compute_shape_function_values (const std::vector< Point< dim >> &unit_points) |
const double & | shape (const unsigned int qpoint, const unsigned int shape_nr) const |
double & | shape (const unsigned int qpoint, const unsigned int shape_nr) |
const Tensor< 1, dim > & | derivative (const unsigned int qpoint, const unsigned int shape_nr) const |
Tensor< 1, dim > & | derivative (const unsigned int qpoint, const unsigned int shape_nr) |
const Tensor< 2, dim > & | second_derivative (const unsigned int qpoint, const unsigned int shape_nr) const |
Tensor< 2, dim > & | second_derivative (const unsigned int qpoint, const unsigned int shape_nr) |
const Tensor< 3, dim > & | third_derivative (const unsigned int qpoint, const unsigned int shape_nr) const |
Tensor< 3, dim > & | third_derivative (const unsigned int qpoint, const unsigned int shape_nr) |
const Tensor< 4, dim > & | fourth_derivative (const unsigned int qpoint, const unsigned int shape_nr) const |
Tensor< 4, dim > & | fourth_derivative (const unsigned int qpoint, const unsigned int shape_nr) |
virtual std::size_t | memory_consumption () const override |
Storage for internal data of polynomial mappings. See Mapping::InternalDataBase for an extensive description.
For the current class, the InternalData class stores data that is computed once when the object is created (in get_data()) as well as data the class wants to store from between the call to fill_fe_values(), fill_fe_face_values(), or fill_fe_subface_values() until possible later calls from the finite element to functions such as transform(). The latter class of member variables are marked as 'mutable'.
Definition at line 274 of file mapping_q.h.
using MappingQ< dim, spacedim >::InternalData::VectorizedArrayType = VectorizedArray<double, std::min<std::size_t>(VectorizedArray<double>::size(), (dim <= 2 ? 2 : 4))> |
For the fast tensor-product path of the MappingQ class, we choose SIMD vectors that are as wide as possible to minimize the number of arithmetic operations. However, we do not want to choose it wider than necessary, e.g., we avoid something like 8-wide AVX-512 when we only compute 3 components of a 3d computation. This is because the additional lanes would not do useful work, but a few operations on very wide vectors can already lead to a lower clock frequency of processors over long time spans (thousands of clock cycles). Hence, we choose 2-wide SIMD for 1D and 2d and 4-wide SIMD for 3d. Note that we do not immediately fall back to no SIMD for 1d because all architectures that support SIMD also support 128-bit vectors (and none is reported to reduce clock frequency for 128-bit SIMD).
Definition at line 490 of file mapping_q.h.
MappingQ< dim, spacedim >::InternalData::InternalData | ( | const unsigned int | polynomial_degree | ) |
Constructor. The argument denotes the polynomial degree of the mapping to which this object will correspond.
Definition at line 51 of file mapping_q.cc.
void MappingQ< dim, spacedim >::InternalData::initialize | ( | const UpdateFlags | update_flags, |
const Quadrature< dim > & | quadrature, | ||
const unsigned int | n_original_q_points | ||
) |
Initialize the object's member variables related to cell data based on the given arguments.
The function also calls compute_shape_function_values() to actually set the member variables related to the values and derivatives of the mapping shape functions.
Definition at line 84 of file mapping_q.cc.
void MappingQ< dim, spacedim >::InternalData::initialize_face | ( | const UpdateFlags | update_flags, |
const Quadrature< dim > & | quadrature, | ||
const unsigned int | n_original_q_points | ||
) |
Initialize the object's member variables related to cell and face data based on the given arguments. In order to initialize cell data, this function calls initialize().
Definition at line 214 of file mapping_q.cc.
void MappingQ< dim, spacedim >::InternalData::compute_shape_function_values | ( | const std::vector< Point< dim >> & | unit_points | ) |
Compute the values and/or derivatives of the shape functions used for the mapping.
Which values, derivatives, or higher order derivatives are computed is determined by which of the member arrays have nonzero sizes. They are typically set to their appropriate sizes by the initialize() and initialize_face() functions, which indeed call this function internally. However, it is possible (and at times useful) to do the resizing by hand and then call this function directly. An example is in a Newton iteration where we update the location of a quadrature point (e.g., in MappingQ::transform_real_to_uni_cell()) and need to re- compute the mapping and its derivatives at this location, but have already sized all internal arrays correctly.
Definition at line 270 of file mapping_q.cc.
const double& MappingQ< dim, spacedim >::InternalData::shape | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) | const |
Shape function at quadrature point. Shape functions are in tensor product order, so vertices must be reordered to obtain transformation.
double& MappingQ< dim, spacedim >::InternalData::shape | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) |
Shape function at quadrature point. See above.
const Tensor<1, dim>& MappingQ< dim, spacedim >::InternalData::derivative | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) | const |
Gradient of shape function in quadrature point. See above.
Tensor<1, dim>& MappingQ< dim, spacedim >::InternalData::derivative | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) |
Gradient of shape function in quadrature point. See above.
const Tensor<2, dim>& MappingQ< dim, spacedim >::InternalData::second_derivative | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) | const |
Second derivative of shape function in quadrature point. See above.
Tensor<2, dim>& MappingQ< dim, spacedim >::InternalData::second_derivative | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) |
Second derivative of shape function in quadrature point. See above.
const Tensor<3, dim>& MappingQ< dim, spacedim >::InternalData::third_derivative | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) | const |
third derivative of shape function in quadrature point. See above.
Tensor<3, dim>& MappingQ< dim, spacedim >::InternalData::third_derivative | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) |
third derivative of shape function in quadrature point. See above.
const Tensor<4, dim>& MappingQ< dim, spacedim >::InternalData::fourth_derivative | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) | const |
fourth derivative of shape function in quadrature point. See above.
Tensor<4, dim>& MappingQ< dim, spacedim >::InternalData::fourth_derivative | ( | const unsigned int | qpoint, |
const unsigned int | shape_nr | ||
) |
fourth derivative of shape function in quadrature point. See above.
|
overridevirtual |
Return an estimate (in bytes) for the memory consumption of this object.
Definition at line 63 of file mapping_q.cc.
AlignedVector<double> MappingQ< dim, spacedim >::InternalData::shape_values |
Values of shape functions. Access by function shape
.
Computed once.
Definition at line 399 of file mapping_q.h.
AlignedVector<Tensor<1, dim> > MappingQ< dim, spacedim >::InternalData::shape_derivatives |
Values of shape function derivatives. Access by function derivative
.
Computed once.
Definition at line 406 of file mapping_q.h.
AlignedVector<Tensor<2, dim> > MappingQ< dim, spacedim >::InternalData::shape_second_derivatives |
Values of shape function second derivatives. Access by function second_derivative
.
Computed once.
Definition at line 414 of file mapping_q.h.
AlignedVector<Tensor<3, dim> > MappingQ< dim, spacedim >::InternalData::shape_third_derivatives |
Values of shape function third derivatives. Access by function second_derivative
.
Computed once.
Definition at line 422 of file mapping_q.h.
AlignedVector<Tensor<4, dim> > MappingQ< dim, spacedim >::InternalData::shape_fourth_derivatives |
Values of shape function fourth derivatives. Access by function second_derivative
.
Computed once.
Definition at line 430 of file mapping_q.h.
std::array<std::vector<Tensor<1, dim> >, GeometryInfo<dim>::faces_per_cell *(dim - 1)> MappingQ< dim, spacedim >::InternalData::unit_tangentials |
Unit tangential vectors. Used for the computation of boundary forms and normal vectors.
This array has (dim-1) * GeometryInfo::faces_per_cell
entries. The first GeometryInfo::faces_per_cell contain the vectors in the first tangential direction for each face; the second set of GeometryInfo::faces_per_cell entries contain the vectors in the second tangential direction (only in 3d, since there we have 2 tangential directions per face), etc.
Filled once.
Definition at line 447 of file mapping_q.h.
const unsigned int MappingQ< dim, spacedim >::InternalData::polynomial_degree |
The polynomial degree of the mapping. Since the objects here are also used (with minor adjustments) by MappingQ, we need to store this.
Definition at line 453 of file mapping_q.h.
const unsigned int MappingQ< dim, spacedim >::InternalData::n_shape_functions |
Number of shape functions. If this is a Q1 mapping, then it is simply the number of vertices per cell. However, since also derived classes use this class (e.g. the Mapping_Q() class), the number of shape functions may also be different.
In general, it is \((p+1)^\text{dim}\), where \(p\) is the polynomial degree of the mapping.
Definition at line 464 of file mapping_q.h.
QGaussLobatto<1> MappingQ< dim, spacedim >::InternalData::line_support_points |
Definition at line 474 of file mapping_q.h.
internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType> MappingQ< dim, spacedim >::InternalData::shape_info |
In case the quadrature rule given represents a tensor product we need to store the evaluations of the 1d polynomials at the 1d quadrature points. That is what this variable is for.
Definition at line 500 of file mapping_q.h.
|
mutable |
In case the quadrature rule given represents a tensor product we need to store temporary data in this object.
Definition at line 506 of file mapping_q.h.
bool MappingQ< dim, spacedim >::InternalData::tensor_product_quadrature |
Indicates whether the given Quadrature object is a tensor product.
Definition at line 511 of file mapping_q.h.
|
mutable |
Tensors of covariant transformation at each of the quadrature points. The matrix stored is the Jacobian * G^{-1}, where G = Jacobian^{t} * Jacobian, is the first fundamental form of the map; if dim=spacedim then it reduces to the transpose of the inverse of the Jacobian matrix, which itself is stored in the contravariant
field of this structure.
Computed on each cell.
Definition at line 522 of file mapping_q.h.
|
mutable |
Tensors of contravariant transformation at each of the quadrature points. The contravariant matrix is the Jacobian of the transformation, i.e. \(J_{ij}=dx_i/d\hat x_j\).
Computed on each cell.
Definition at line 531 of file mapping_q.h.
|
mutable |
Auxiliary vectors for internal use.
Definition at line 536 of file mapping_q.h.
|
mutable |
Stores the support points of the mapping shape functions on the cell_of_current_support_points
.
Definition at line 542 of file mapping_q.h.
|
mutable |
Stores the cell of which the mapping_support_points
are stored.
Definition at line 548 of file mapping_q.h.
|
mutable |
The determinant of the Jacobian in each quadrature point. Filled if update_volume_elements.
Definition at line 554 of file mapping_q.h.
|
inherited |
A set of update flags specifying the kind of information that an implementation of the Mapping interface needs to compute on each cell or face, i.e., in Mapping::fill_fe_values() and friends.
This set of flags is stored here by implementations of Mapping::get_data(), Mapping::get_face_data(), or Mapping::get_subface_data(), and is that subset of the update flags passed to those functions that require re-computation on every cell. (The subset of the flags corresponding to information that can be computed once and for all already at the time of the call to Mapping::get_data() – or an implementation of that interface – need not be stored here because it has already been taken care of.)