deal.II version GIT relicensing-2165-gc91f007519 2024-11-20 01:40:00+00:00
|
#include <deal.II/matrix_free/mapping_info.h>
Public Member Functions | |
void | initialize (const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const FaceInfo< VectorizedArrayType::size()> &face_info, const std::vector< unsigned int > &active_fe_index, const std::shared_ptr<::hp::MappingCollection< dim > > &mapping, const std::vector<::hp::QCollection< dim > > &quad, const UpdateFlags update_flags_cells, const UpdateFlags update_flags_boundary_faces, const UpdateFlags update_flags_inner_faces, const UpdateFlags update_flags_faces_by_cells, const bool piola_transform) |
void | update_mapping (const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const FaceInfo< VectorizedArrayType::size()> &face_info, const std::vector< unsigned int > &active_fe_index, const std::shared_ptr<::hp::MappingCollection< dim > > &mapping) |
GeometryType | get_cell_type (const unsigned int cell_chunk_no) const |
void | clear () |
std::size_t | memory_consumption () const |
template<typename StreamType > | |
void | print_memory_consumption (StreamType &out, const TaskInfo &task_info) const |
void | compute_mapping_q (const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const FaceInfo< VectorizedArrayType::size()> &face_info) |
void | initialize_cells (const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const std::vector< unsigned int > &active_fe_index, const ::hp::MappingCollection< dim > &mapping) |
void | initialize_faces (const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const std::vector< FaceToCellTopology< VectorizedArrayType::size()> > &faces, const std::vector< unsigned int > &active_fe_index, const ::hp::MappingCollection< dim > &mapping) |
void | initialize_faces_by_cells (const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const FaceInfo< VectorizedArrayType::size()> &face_info, const ::hp::MappingCollection< dim > &mapping) |
Public Attributes | |
UpdateFlags | update_flags_cells |
UpdateFlags | update_flags_boundary_faces |
UpdateFlags | update_flags_inner_faces |
UpdateFlags | update_flags_faces_by_cells |
std::vector< GeometryType > | cell_type |
std::vector< GeometryType > | face_type |
std::vector< std::array< GeometryType, GeometryInfo< dim >::faces_per_cell > > | faces_by_cells_type |
std::vector< MappingInfoStorage< dim, dim, VectorizedArrayType > > | cell_data |
std::vector< MappingInfoStorage< dim - 1, dim, VectorizedArrayType > > | face_data |
std::vector< MappingInfoStorage< dim - 1, dim, VectorizedArrayType > > | face_data_by_cells |
std::shared_ptr<::hp::MappingCollection< dim > > | mapping_collection |
ObserverPointer< const Mapping< dim > > | mapping |
std::vector< std::vector< ReferenceCell > > | reference_cell_types |
The class that stores all geometry-dependent data related with cell interiors for use in the matrix-free class.
Definition at line 54 of file mapping_info.h.
void internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::initialize | ( | const ::Triangulation< dim > & | tria, |
const std::vector< std::pair< unsigned int, unsigned int > > & | cells, | ||
const FaceInfo< VectorizedArrayType::size()> & | face_info, | ||
const std::vector< unsigned int > & | active_fe_index, | ||
const std::shared_ptr<::hp::MappingCollection< dim > > & | mapping, | ||
const std::vector<::hp::QCollection< dim > > & | quad, | ||
const UpdateFlags | update_flags_cells, | ||
const UpdateFlags | update_flags_boundary_faces, | ||
const UpdateFlags | update_flags_inner_faces, | ||
const UpdateFlags | update_flags_faces_by_cells, | ||
const bool | piola_transform | ||
) |
Compute the information in the given cells and faces. The cells are specified by the level and the index within the level (as given by CellIterator::level() and CellIterator::index(), in order to allow for different kinds of iterators, e.g. standard DoFHandler, multigrid, etc.) on a fixed Triangulation. In addition, a mapping and several 1d quadrature formulas are given.
void internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::update_mapping | ( | const ::Triangulation< dim > & | tria, |
const std::vector< std::pair< unsigned int, unsigned int > > & | cells, | ||
const FaceInfo< VectorizedArrayType::size()> & | face_info, | ||
const std::vector< unsigned int > & | active_fe_index, | ||
const std::shared_ptr<::hp::MappingCollection< dim > > & | mapping | ||
) |
Update the information in the given cells and faces that is the result of a change in the given mapping
class, keeping the cells, quadrature formulas and other unknowns unchanged. This call is only valid if MappingInfo::initialize() has been called before.
|
inline |
Return the type of a given cell as detected during initialization.
Definition at line 304 of file mapping_info.h.
void internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::clear | ( | ) |
Clear all data fields in this class.
std::size_t internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::memory_consumption | ( | ) | const |
Return the memory consumption of this class in bytes.
void internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::print_memory_consumption | ( | StreamType & | out, |
const TaskInfo & | task_info | ||
) | const |
Prints a detailed summary of memory consumption in the different structures of this class to the given output stream.
void internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::compute_mapping_q | ( | const ::Triangulation< dim > & | tria, |
const std::vector< std::pair< unsigned int, unsigned int > > & | cells, | ||
const FaceInfo< VectorizedArrayType::size()> & | face_info | ||
) |
Internal function to compute the geometry for the case the mapping is a MappingQ and a single quadrature formula per slot (non-hp-case) is used. This method computes all data from the underlying cell quadrature points using the fast operator evaluation techniques from the matrix-free framework itself, i.e., it uses a polynomial description of the cell geometry (that is computed in a first step) and then computes all Jacobians and normal vectors based on this information. This optimized approach is much faster than going through FEValues and FEFaceValues, especially when several different quadrature formulas are involved, and consumes less memory.
tria | The triangulation to be used for setup |
cells | The actual cells of the triangulation to be worked on, given as a tuple of the level and index within the level as used in the main initialization of the class |
face_info | The description of the connectivity from faces to cells as filled in the MatrixFree class |
void internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::initialize_cells | ( | const ::Triangulation< dim > & | tria, |
const std::vector< std::pair< unsigned int, unsigned int > > & | cells, | ||
const std::vector< unsigned int > & | active_fe_index, | ||
const ::hp::MappingCollection< dim > & | mapping | ||
) |
Computes the information in the given cells, called within initialize.
void internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::initialize_faces | ( | const ::Triangulation< dim > & | tria, |
const std::vector< std::pair< unsigned int, unsigned int > > & | cells, | ||
const std::vector< FaceToCellTopology< VectorizedArrayType::size()> > & | faces, | ||
const std::vector< unsigned int > & | active_fe_index, | ||
const ::hp::MappingCollection< dim > & | mapping | ||
) |
Computes the information in the given faces, called within initialize.
void internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::initialize_faces_by_cells | ( | const ::Triangulation< dim > & | tria, |
const std::vector< std::pair< unsigned int, unsigned int > > & | cells, | ||
const FaceInfo< VectorizedArrayType::size()> & | face_info, | ||
const ::hp::MappingCollection< dim > & | mapping | ||
) |
Computes the information in the given faces, called within initialize.
UpdateFlags internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::update_flags_cells |
The given update flags for computing the geometry on the cells.
Definition at line 122 of file mapping_info.h.
UpdateFlags internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::update_flags_boundary_faces |
The given update flags for computing the geometry on the boundary faces.
Definition at line 128 of file mapping_info.h.
UpdateFlags internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::update_flags_inner_faces |
The given update flags for computing the geometry on the interior faces.
Definition at line 134 of file mapping_info.h.
UpdateFlags internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::update_flags_faces_by_cells |
The given update flags for computing the geometry on the faces for cell-centric loops.
Definition at line 140 of file mapping_info.h.
std::vector<GeometryType> internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::cell_type |
Stores whether a cell is Cartesian (cell type 0), has constant transform data (Jacobians) (cell type 1), or is general (cell type 3). Type 2 is only used for faces and no cells are assigned this value.
Definition at line 148 of file mapping_info.h.
std::vector<GeometryType> internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::face_type |
Stores whether a face (and both cells adjacent to the face) is Cartesian (face type 0), whether it represents an affine situation (face type 1), whether it is a flat face where the normal vector is the same throughout the face (face type 2), or is general (face type 3).
Definition at line 157 of file mapping_info.h.
std::vector<std::array<GeometryType, GeometryInfo<dim>::faces_per_cell> > internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::faces_by_cells_type |
Store whether the face geometry for the face_data_by_cells data structure is represented as Cartesian (cell type 0), with constant transform data (Jacobians) (cell type 1), or a general type (cell type 3). Note that both the interior and exterior agree on the type of the data structure, using the more general of the two.
Definition at line 167 of file mapping_info.h.
std::vector<MappingInfoStorage<dim, dim, VectorizedArrayType> > internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::cell_data |
The data cache for the cells.
Definition at line 172 of file mapping_info.h.
std::vector<MappingInfoStorage<dim - 1, dim, VectorizedArrayType> > internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::face_data |
The data cache for the faces.
Definition at line 178 of file mapping_info.h.
std::vector<MappingInfoStorage<dim - 1, dim, VectorizedArrayType> > internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::face_data_by_cells |
The data cache for the face-associated-with-cell topology, following the cell_type
variable for the cell types.
Definition at line 185 of file mapping_info.h.
std::shared_ptr<::hp::MappingCollection<dim> > internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::mapping_collection |
The pointer to the underlying hp::MappingCollection object.
Definition at line 190 of file mapping_info.h.
ObserverPointer<const Mapping<dim> > internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::mapping |
The pointer to the first entry of mapping_collection.
Definition at line 195 of file mapping_info.h.
std::vector<std::vector<ReferenceCell> > internal::MatrixFreeFunctions::MappingInfo< dim, Number, VectorizedArrayType >::reference_cell_types |
Reference-cell type related to each quadrature and active quadrature index.
Definition at line 201 of file mapping_info.h.