17#ifndef dealii_matrix_free_mapping_info_h
18#define dealii_matrix_free_mapping_info_h
46 namespace MatrixFreeFunctions
114 template <
int structdim,
117 typename VectorizedArrayType>
121 std::is_same<Number, typename VectorizedArrayType::value_type>::value,
122 "Type of Number and of VectorizedArrayType do not match.");
242 std::array<AlignedVector<Tensor<2, spacedim, VectorizedArrayType>>, 2>
262 spacedim *(spacedim + 1) / 2,
274 std::array<AlignedVector<Tensor<1, spacedim, VectorizedArrayType>>, 2>
313 template <
typename StreamType>
333 template <
int dim,
typename Number,
typename VectorizedArrayType>
346 const ::Triangulation<dim> & tria,
347 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
348 const FaceInfo<VectorizedArrayType::size()> & faces,
349 const std::vector<unsigned int> &active_fe_index,
365 const ::Triangulation<dim> & tria,
366 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
367 const FaceInfo<VectorizedArrayType::size()> & faces,
368 const std::vector<unsigned int> &active_fe_index,
393 template <
typename StreamType>
441 std::vector<MappingInfoStorage<dim, dim, Number, VectorizedArrayType>>
496 const ::Triangulation<dim> & tria,
497 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
507 const ::Triangulation<dim> & tria,
508 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
509 const std::vector<unsigned int> & active_fe_index,
510 const ::hp::MappingCollection<dim> &
mapping);
518 const ::Triangulation<dim> & tria,
519 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
522 const std::vector<unsigned int> & active_fe_index,
523 const ::hp::MappingCollection<dim> &
mapping);
531 const ::Triangulation<dim> & tria,
532 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
533 const ::hp::MappingCollection<dim> &
mapping);
552 template <
int,
typename,
bool,
typename>
555 template <
int dim,
typename Number,
typename VectorizedArrayType>
560 const unsigned int quad_no)
563 return &mapping_info.cell_data[quad_no];
567 template <
int dim,
typename Number,
typename VectorizedArrayType>
573 const unsigned int quad_no)
576 return &mapping_info.face_data[quad_no];
593 template <
typename Number,
604 const std::vector<Number> &v2)
const;
612 const Tensor<1, VectorizedArrayType::size(), Number> &t1,
613 const Tensor<1, VectorizedArrayType::size(), Number> &t2)
const;
622 const Tensor<1, dim,
Tensor<1, VectorizedArrayType::size(), Number>>
624 const Tensor<1, dim,
Tensor<1, VectorizedArrayType::size(), Number>>
634 const Tensor<2, dim,
Tensor<1, VectorizedArrayType::size(), Number>>
636 const Tensor<2, dim,
Tensor<1, VectorizedArrayType::size(), Number>>
654 template <
int structdim,
657 typename VectorizedArrayType>
662 for (
unsigned int i = 0; i < descriptor.size(); ++i)
663 if (n_q_points == descriptor[i].n_q_points)
670 template <
int dim,
typename Number,
typename VectorizedArrayType>
673 const unsigned int cell_no)
const
676 return cell_type[cell_no];
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
@ update_default
No update.
#define AssertIndexRange(index, range)
bool operator()(const std::vector< Number > &v1, const std::vector< Number > &v2) const
FPArrayComparator(const Number scaling)
bool operator()(const std::array< Tensor< 2, dim, Number >, dim+1 > &t1, const std::array< Tensor< 2, dim, Number >, dim+1 > &t2) const
bool operator()(const Tensor< 1, dim, Tensor< 1, VectorizedArrayType::size(), Number > > &t1, const Tensor< 1, dim, Tensor< 1, VectorizedArrayType::size(), Number > > &t2) const
bool operator()(const Tensor< 1, VectorizedArrayType::size(), Number > &t1, const Tensor< 1, VectorizedArrayType::size(), Number > &t2) const
bool operator()(const Tensor< 2, dim, Tensor< 1, VectorizedArrayType::size(), Number > > &t1, const Tensor< 2, dim, Tensor< 1, VectorizedArrayType::size(), Number > > &t2) const
static const MappingInfoStorage< dim, dim, Number, VectorizedArrayType > * get(const MappingInfo< dim, Number, VectorizedArrayType > &mapping_info, const unsigned int quad_no)
static const MappingInfoStorage< dim - 1, dim, Number, VectorizedArrayType > * get(const MappingInfo< dim, Number, VectorizedArrayType > &mapping_info, const unsigned int quad_no)
Quadrature< structdim > quadrature
Quadrature< 1 > quadrature_1d
std::size_t memory_consumption() const
AlignedVector< Number > quadrature_weights
void initialize(const Quadrature< dim_q > &quadrature, const UpdateFlags update_flags_inner_faces=update_default)
::Table< 2, unsigned int > face_orientations
void initialize(const Quadrature< 1 > &quadrature_1d, const UpdateFlags update_flags_inner_faces=update_default)
std::array< AlignedVector< Number >, structdim > tensor_quadrature_weights
std::array< AlignedVector< Tensor< 1, spacedim *(spacedim+1)/2, Tensor< 1, spacedim, VectorizedArrayType > > >, 2 > jacobian_gradients
AlignedVector< Point< spacedim, VectorizedArrayType > > quadrature_points
unsigned int quad_index_from_n_q_points(const unsigned int n_q_points) const
std::vector<::hp::QCollection< structdim > > q_collection
AlignedVector< unsigned int > quadrature_point_offsets
AlignedVector< Tensor< 1, spacedim, VectorizedArrayType > > normal_vectors
std::size_t memory_consumption() const
AlignedVector< VectorizedArrayType > JxW_values
AlignedVector< unsigned int > data_index_offsets
std::array< AlignedVector< Tensor< 1, spacedim, VectorizedArrayType > >, 2 > normals_times_jacobians
void print_memory_consumption(StreamType &out, const TaskInfo &task_info) const
std::array< AlignedVector< Tensor< 2, spacedim, VectorizedArrayType > >, 2 > jacobians
std::vector< QuadratureDescriptor > descriptor
std::vector< MappingInfoStorage< dim, dim, Number, VectorizedArrayType > > cell_data
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)
std::vector< MappingInfoStorage< dim - 1, dim, Number, VectorizedArrayType > > face_data_by_cells
std::shared_ptr<::hp::MappingCollection< dim > > mapping_collection
GeometryType get_cell_type(const unsigned int cell_chunk_no) const
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)
std::vector< GeometryType > cell_type
UpdateFlags update_flags_inner_faces
void initialize_faces_by_cells(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const ::hp::MappingCollection< dim > &mapping)
std::vector< GeometryType > face_type
std::size_t memory_consumption() const
UpdateFlags update_flags_boundary_faces
std::vector< MappingInfoStorage< dim - 1, dim, Number, VectorizedArrayType > > face_data
std::vector< std::vector<::ReferenceCell > > reference_cell_types
UpdateFlags update_flags_cells
UpdateFlags update_flags_faces_by_cells
SmartPointer< const Mapping< dim > > mapping
void initialize(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const FaceInfo< VectorizedArrayType::size()> &faces, 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)
void update_mapping(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const FaceInfo< VectorizedArrayType::size()> &faces, const std::vector< unsigned int > &active_fe_index, const std::shared_ptr<::hp::MappingCollection< dim > > &mapping)
void print_memory_consumption(StreamType &out, const TaskInfo &task_info) const
static UpdateFlags compute_update_flags(const UpdateFlags update_flags, const std::vector<::hp::QCollection< dim > > &quad=std::vector<::hp::QCollection< dim > >())
void compute_mapping_q(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const std::vector< FaceToCellTopology< VectorizedArrayType::size()> > &faces)