17 #ifndef dealii_matrix_free_fe_evaluation_data_h
18 #define dealii_matrix_free_fe_evaluation_data_h
49 <<
"You are requesting information from an FEEvaluation/FEFaceEvaluation "
50 <<
"object for which this kind of information has not been computed. What "
51 <<
"information these objects compute is determined by the update_* flags you "
52 <<
"pass to MatrixFree::reinit() via MatrixFree::AdditionalData. Here, "
53 <<
"the operation you are attempting requires the <" << arg1
54 <<
"> flag to be set, but it was apparently not specified "
55 <<
"upon initialization.");
61 int n_q_points_1d = fe_degree + 1,
62 int n_components_ = 1,
76 namespace MatrixFreeFunctions
78 template <
int,
typename>
79 class MappingDataOnTheFly;
112 template <
int dim,
typename Number,
bool is_face>
117 MappingInfoStorage<(is_face ? dim - 1 : dim), dim, Number>;
158 const unsigned int n_components);
187 JxW(
const unsigned int q_point)
const;
381 const std::vector<unsigned int> &
479 const std::array<unsigned int, n_lanes> &
493 const std::array<unsigned int, n_lanes> &
522 const std::array<unsigned int, n_lanes> &
574 template <
typename T>
575 std::array<
T, Number::size()>
583 template <
typename T>
586 const std::array<
T, Number::size()> & value)
const;
614 template <
typename T>
615 std::array<
T, Number::size()>
623 template <
typename T>
626 const std::array<
T, Number::size()> & value)
const;
661 const std::shared_ptr<
985 template <
int,
int,
typename,
bool,
typename>
988 template <
int,
int,
int,
int,
typename,
typename>
1009 template <
int dim,
typename Number,
bool is_face>
1014 InitializationData{&shape_info, nullptr, nullptr, 0, 0, nullptr},
1021 template <
int dim,
typename Number,
bool is_face>
1023 const InitializationData &initialization_data,
1027 :
data(initialization_data.shape_info)
1070 template <
int dim,
typename Number,
bool is_face>
1072 const std::shared_ptr<
1113 template <
int dim,
typename Number,
bool is_face>
1164 template <
int dim,
typename Number,
bool is_face>
1172 const unsigned int tensor_dofs_per_component =
1173 Utilities::fixed_power<dim>(
data->
data.front().fe_degree + 1);
1176 const unsigned int size_scratch_data =
1180 const unsigned int size_data_arrays =
1182 (
n_components * ((dim * (dim + 1)) / 2 + 2 * dim + 1) *
1185 const unsigned int allocated_size = size_scratch_data + size_data_arrays;
1205 template <
int dim,
typename Number,
bool is_face>
1211 ExcMessage(
"Faces can only be set if the is_face template parameter "
1236 template <
int dim,
typename Number,
bool is_face>
1239 const unsigned int q_point)
const
1244 "update_normal_vectors"));
1253 template <
int dim,
typename Number,
bool is_face>
1260 "update_values|update_gradients"));
1272 template <
int dim,
typename Number,
bool is_face>
1275 const unsigned int q)
const
1280 "update_quadrature_points"));
1285 if (is_face ==
false &&
1293 for (
unsigned int d = 0;
d < dim; ++
d)
1294 point[
d] += jac[
d][
d] *
static_cast<typename Number::value_type
>(
1297 for (
unsigned int d = 0;
d < dim; ++
d)
1298 for (
unsigned int e = 0;
e < dim; ++
e)
1299 point[
d] += jac[
d][
e] *
static_cast<typename Number::value_type
>(
1309 template <
int dim,
typename Number,
bool is_face>
1312 const unsigned int q_point)
const
1317 "update_gradients"));
1326 template <
int dim,
typename Number,
bool is_face>
1327 inline const Number *
1335 template <
int dim,
typename Number,
bool is_face>
1347 template <
int dim,
typename Number,
bool is_face>
1348 inline const Number *
1359 template <
int dim,
typename Number,
bool is_face>
1372 template <
int dim,
typename Number,
bool is_face>
1373 inline const Number *
1385 template <
int dim,
typename Number,
bool is_face>
1398 template <
int dim,
typename Number,
bool is_face>
1399 inline const Number *
1410 template <
int dim,
typename Number,
bool is_face>
1422 template <
int dim,
typename Number,
bool is_face>
1439 template <
int dim,
typename Number,
bool is_face>
1451 template <
int dim,
typename Number,
bool is_face>
1461 template <
int dim,
typename Number,
bool is_face>
1467 "FEEvaluation was not initialized with a MatrixFree object!"));
1473 template <
int dim,
typename Number,
bool is_face>
1474 inline const std::vector<unsigned int> &
1482 template <
int dim,
typename Number,
bool is_face>
1491 template <
int dim,
typename Number,
bool is_face>
1504 template <
int dim,
typename Number,
bool is_face>
1513 template <
int dim,
typename Number,
bool is_face>
1522 template <
int dim,
typename Number,
bool is_face>
1531 template <
int dim,
typename Number,
bool is_face>
1540 template <
int dim,
typename Number,
bool is_face>
1549 ExcMessage(
"All face numbers can only be queried for ECL at exterior "
1550 "faces. Use get_face_no() in other cases."));
1557 template <
int dim,
typename Number,
bool is_face>
1566 template <
int dim,
typename Number,
bool is_face>
1569 const unsigned int v)
const
1576 ExcMessage(
"All face numbers can only be queried for ECL at exterior "
1577 "faces. Use get_face_no() in other cases."));
1584 template <
int dim,
typename Number,
bool is_face>
1593 template <
int dim,
typename Number,
bool is_face>
1602 template <
int dim,
typename Number,
bool is_face>
1613 template <std::size_t
N,
1614 typename VectorizedArrayType2,
1615 typename GlobalVectorType,
1618 process_cell_or_face_data(
const std::array<unsigned int, N> indices,
1619 GlobalVectorType & array,
1620 VectorizedArrayType2 & out,
1623 for (
unsigned int i = 0; i <
N; ++i)
1627 fu(out[i], array[indices[i] /
N][indices[i] %
N]);
1634 template <
int dim,
typename Number,
bool is_face>
1639 Number out = Number(1.);
1640 internal::process_cell_or_face_data(this->get_cell_ids(),
1643 [](
auto &local,
const auto &global) {
1651 template <
int dim,
typename Number,
bool is_face>
1655 const Number & in)
const
1657 internal::process_cell_or_face_data(this->get_cell_ids(),
1660 [](
const auto &local,
auto &global) {
1667 template <
int dim,
typename Number,
bool is_face>
1668 template <
typename T>
1669 inline std::array<
T, Number::size()>
1673 std::array<
T, Number::size()> out;
1674 internal::process_cell_or_face_data(this->get_cell_ids(),
1677 [](
auto &local,
const auto &global) {
1685 template <
int dim,
typename Number,
bool is_face>
1686 template <
typename T>
1690 const std::array<
T, Number::size()> & in)
const
1692 internal::process_cell_or_face_data(this->get_cell_ids(),
1695 [](
const auto &local,
auto &global) {
1702 template <
int dim,
typename Number,
bool is_face>
1707 Number out = Number(1.);
1708 internal::process_cell_or_face_data(this->get_face_ids(),
1711 [](
auto &local,
const auto &global) {
1719 template <
int dim,
typename Number,
bool is_face>
1723 const Number & in)
const
1725 internal::process_cell_or_face_data(this->get_face_ids(),
1728 [](
const auto &local,
auto &global) {
1735 template <
int dim,
typename Number,
bool is_face>
1736 template <
typename T>
1737 inline std::array<
T, Number::size()>
1741 std::array<
T, Number::size()> out;
1742 internal::process_cell_or_face_data(this->get_face_ids(),
1745 [](
auto &local,
const auto &global) {
1753 template <
int dim,
typename Number,
bool is_face>
1754 template <
typename T>
1758 const std::array<
T, Number::size()> & in)
const
1760 internal::process_cell_or_face_data(this->get_face_ids(),
1763 [](
const auto &local,
auto &global) {
void resize_fast(const size_type new_size)
void reinit(value_type *starting_element, const std::size_t n_elements)
AlignedVector< VectorizedArrayType > * scratch_data_array
const unsigned int quad_no
bool hessians_quad_submitted
internal::MatrixFreeFunctions::GeometryType get_cell_type() const
const Tensor< 2, dim, Number > * jacobian
const MappingInfoStorageType::QuadratureDescriptor * descriptor
Point< dim, Number > quadrature_point(const unsigned int q) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
const unsigned int n_quadrature_points
const MappingInfoStorageType * mapping_data
internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex dof_access_index
std::uint8_t get_face_no(const unsigned int v=0) const
ArrayView< Number > scratch_data
internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex get_dof_access_index() const
const Point< dim, Number > * quadrature_points
unsigned int subface_index
bool values_quad_submitted
FEEvaluationData(const FEEvaluationData &other)=default
void reinit_face(const internal::MatrixFreeFunctions::FaceToCellTopology< n_lanes > &face)
Number read_face_data(const AlignedVector< Number > &array) const
Number JxW(const unsigned int q_point) const
FEEvaluationData & operator=(const FEEvaluationData &other)
const std::array< unsigned int, n_lanes > & get_cell_or_face_ids() const
unsigned int get_first_selected_component() const
const unsigned int n_fe_components
void set_cell_data(AlignedVector< Number > &array, const Number &value) const
const ShapeInfoType * data
Number * gradients_from_hessians_quad
std::shared_ptr< internal::MatrixFreeFunctions::MappingDataOnTheFly< dim, Number > > mapped_geometry
internal::MatrixFreeFunctions::DoFInfo DoFInfo
unsigned int get_active_quadrature_index() const
unsigned int get_mapping_data_index_offset() const
bool gradients_quad_initialized
void set_data_pointers(AlignedVector< Number > *scratch_data, const unsigned int n_components)
std::array< std::uint8_t, n_lanes > face_numbers
const unsigned int active_fe_index
static constexpr unsigned int n_lanes
const Tensor< 1, dim, Number > * normal_vectors
Tensor< 1, dim, Number > get_normal_vector(const unsigned int q_point) const
const Number * begin_values() const
internal::MatrixFreeFunctions::GeometryType cell_type
std::array< T, Number::size()> read_face_data(const AlignedVector< std::array< T, Number::size()>> &array) const
const Number * begin_hessians() const
const std::array< unsigned int, n_lanes > & get_cell_ids() const
Number * begin_dof_values()
std::array< unsigned int, n_lanes > cell_ids
unsigned int get_current_cell_index() const
unsigned int get_subface_index() const
Number * begin_gradients()
Number read_cell_data(const AlignedVector< Number > &array) const
const unsigned int active_quad_index
unsigned int get_active_fe_index() const
const Number * begin_gradients() const
bool gradients_quad_submitted
const std::vector< unsigned int > & get_internal_dof_numbering() const
unsigned int get_quadrature_index() const
std::array< std::uint8_t, n_lanes > face_orientations
const ScalarNumber * quadrature_weights
const Tensor< 1, dim *(dim+1)/2, Tensor< 1, dim, Number > > * jacobian_gradients
std::array< T, Number::size()> read_cell_data(const AlignedVector< std::array< T, Number::size()>> &array) const
const Tensor< 1, dim, Number > * normal_x_jacobian
const internal::MatrixFreeFunctions::DoFInfo & get_dof_info() const
bool is_interior_face() const
const std::array< unsigned int, n_lanes > & get_face_ids() const
ArrayView< Number > get_scratch_data() const
static constexpr unsigned int dimension
internal::MatrixFreeFunctions::ShapeInfo< VectorizedArrayType > ShapeInfoType
FEEvaluationData(const std::shared_ptr< internal::MatrixFreeFunctions::MappingDataOnTheFly< dim, Number >> &mapping_data, const unsigned int n_fe_components, const unsigned int first_selected_component)
Tensor< 2, dim, Number > inverse_jacobian(const unsigned int q_point) const
bool values_quad_initialized
FEEvaluationData(const ShapeInfoType &shape_info, const bool is_interior_face=true)
void set_cell_data(AlignedVector< std::array< T, Number::size()>> &array, const std::array< T, Number::size()> &value) const
bool dof_values_initialized
typename internal::VectorizedArrayTrait< Number >::value_type ScalarNumber
std::array< unsigned int, n_lanes > face_ids
const ShapeInfoType & get_shape_info() const
FEEvaluationData(const InitializationData &initialization_data, const bool is_interior_face, const unsigned int quad_no, const unsigned int first_selected_component)
void set_face_data(AlignedVector< Number > &array, const Number &value) const
Number * begin_hessians()
std::uint8_t get_face_orientation(const unsigned int v=0) const
bool hessians_quad_initialized
unsigned int get_cell_or_face_batch_id() const
void set_face_data(AlignedVector< std::array< T, Number::size()>> &array, const std::array< T, Number::size()> &value) const
const unsigned int first_selected_component
const Number * begin_dof_values() const
const unsigned int dofs_per_component
const unsigned int n_q_points
static constexpr unsigned int n_components
const Point< dim > & point(const unsigned int i) const
#define DEAL_II_ALWAYS_INLINE
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define DeclException0(Exception0)
static ::ExceptionBase & ExcAccessToUninitializedField()
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcNotInitialized()
#define Assert(cond, exc)
static ::ExceptionBase & ExcMatrixFreeAccessToUninitializedMappingField(std::string arg1)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcMessage(std::string arg1)
@ general
No special properties.
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
void quadrature_points(const Triangulation< dim, spacedim > &triangulation, const Quadrature< dim > &quadrature, const std::vector< std::vector< BoundingBox< spacedim >>> &global_bounding_boxes, ParticleHandler< dim, spacedim > &particle_handler, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()), const std::vector< std::vector< double >> &properties={})
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
VectorType::value_type * begin(VectorType &V)
static const unsigned int invalid_unsigned_int
boost::integer_range< IncrementableType > iota_view
const MappingInfoStorageType * mapping_data
unsigned int active_quad_index
const MappingInfoStorageType::QuadratureDescriptor * descriptor
const ShapeInfoType * shape_info
unsigned int active_fe_index
@ dof_access_face_exterior
@ dof_access_face_interior
unsigned char subface_index
unsigned char interior_face_no
std::array< unsigned int, vectorization_width > cells_interior
types::boundary_id exterior_face_no
std::array< unsigned int, vectorization_width > cells_exterior
unsigned char face_orientation
Quadrature< structdim > quadrature
AlignedVector< unsigned int > data_index_offsets
unsigned int dofs_per_component_on_cell
std::vector< UnivariateShapeData< Number > > data
std::vector< unsigned int > lexicographic_numbering