18 #ifndef __deal2__matrix_free_mapping_info_h 19 #define __deal2__matrix_free_mapping_info_h 23 #include <deal.II/base/vectorization.h> 24 #include <deal.II/hp/q_collection.h> 25 #include <deal.II/fe/fe.h> 26 #include <deal.II/fe/mapping.h> 27 #include <deal.II/matrix_free/helper_functions.h> 37 namespace MatrixFreeFunctions
45 template <
int dim,
typename Number>
73 void initialize (const ::Triangulation<dim> &tria,
74 const std::vector<std::pair<unsigned int,unsigned int> > &cells,
75 const std::vector<unsigned int> &active_fe_index,
91 CellType
get_cell_type (
const unsigned int cell_chunk_no)
const;
112 template <
typename STREAM>
270 template <
typename STREAM>
307 void resize (
const unsigned int size);
313 const double jac_size;
320 const std::pair<unsigned int,unsigned int> *cells,
321 const unsigned int cell,
322 const unsigned int my_q,
324 CellType (&cell_t)[VectorizedArray<Number>::n_array_elements],
333 template <
int dim,
typename Number>
347 template <
int dim,
typename Number>
355 return enum_cell_type;
360 template <
int dim,
typename Number>
372 DEAL_II_NAMESPACE_CLOSE
AlignedVector< std::pair< Tensor< 1, dim, VectorizedArray< Number > >, VectorizedArray< Number > > > cartesian_data
::hp::QCollection< dim > quadrature
AlignedVector< VectorizedArray< Number > > JxW_values
AlignedVector< std::pair< Tensor< 2, dim, VectorizedArray< Number > >, VectorizedArray< Number > > > affine_data
bool quadrature_points_initialized
std::vector< unsigned int > n_q_points
#define AssertIndexRange(index, range)
std::vector< MappingInfoDependent > mapping_data_gen
std::vector< unsigned int > rowstart_jacobians
std::vector< unsigned int > n_q_points_face
CellType get_cell_type(const unsigned int cell_chunk_no) const
void print_memory_consumption(STREAM &out, const SizeInfo &size_info) const
UpdateFlags compute_update_flags(const UpdateFlags update_flags, const std::vector<::hp::QCollection< 1 > > &quad) const
std::vector< unsigned int > cell_type
#define Assert(cond, exc)
std::vector< unsigned int > quad_index_conversion
AlignedVector< Tensor< 1,(dim >1?dim *(dim-1)/2:1), Tensor< 1, dim, VectorizedArray< Number > > > > jacobians_grad_upper
void print_memory_consumption(STREAM &out, const SizeInfo &size_info) const
unsigned int quad_index_from_n_q_points(const unsigned int n_q_points) const
std::vector< unsigned int > rowstart_q_points
std::size_t memory_consumption() const
bool JxW_values_initialized
unsigned int get_cell_data_index(const unsigned int cell_chunk_no) const
std::size_t memory_consumption() const
static const unsigned int n_cell_types
::hp::QCollection< dim-1 > face_quadrature
AlignedVector< Point< dim, VectorizedArray< Number > > > quadrature_points
void initialize(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const std::vector< unsigned int > &active_fe_index, const Mapping< dim > &mapping, const std::vector<::hp::QCollection< 1 > > &quad, const UpdateFlags update_flags)
static const std::size_t n_cell_type_bits
bool second_derivatives_initialized
::ExceptionBase & ExcInternalError()
std::vector< AlignedVector< VectorizedArray< Number > > > quadrature_weights
void evaluate_on_cell(const ::Triangulation< dim > &tria, const std::pair< unsigned int, unsigned int > *cells, const unsigned int cell, const unsigned int my_q, CellType(&cell_t_prev)[VectorizedArray< Number >::n_array_elements], CellType(&cell_t)[VectorizedArray< Number >::n_array_elements], FEValues< dim, dim > &fe_values, CellData &cell_data) const
AlignedVector< Tensor< 2, dim, VectorizedArray< Number > > > jacobians
AlignedVector< Tensor< 2, dim, VectorizedArray< Number > > > jacobians_grad_diag