17 #ifndef dealii_non_matching_mapping_info_h
18 #define dealii_non_matching_mapping_info_h
46 template <
int dim,
int spacedim = dim>
59 return mapping->requires_update_flags(update_flags);
70 internal_mapping_data,
82 *internal_mapping_data);
89 internal_mapping_data =
90 mapping->get_data(update_flags_mapping, quadrature);
93 cell_similarity = mapping->fill_fe_values(cell,
96 *internal_mapping_data,
109 internal_mapping_data,
119 *internal_mapping_data);
126 internal_mapping_data =
127 mapping->get_data(update_flags_mapping, quadrature);
130 mapping->fill_fe_immersed_surface_values(cell,
132 *internal_mapping_data,
143 const unsigned int face_no,
146 internal_mapping_data,
157 *internal_mapping_data);
160 ReferenceCells::get_hypercube<dim>(),
163 cell->face_orientation(face_no),
164 cell->face_flip(face_no),
165 cell->face_rotation(face_no)),
168 mapping_q->fill_mapping_data_for_face_quadrature(
169 cell, face_no, quadrature, *internal_mapping_data, mapping_data);
173 auto internal_mapping_data =
174 mapping->get_face_data(update_flags_mapping,
177 mapping->fill_fe_face_values(cell,
180 *internal_mapping_data,
186 template <
int dim,
int spacedim = dim>
193 const auto jac_0 = inverse_jacobians[0];
194 const double zero_tolerance_double =
196 bool jacobian_constant =
true;
197 for (
unsigned int q = 1; q < inverse_jacobians.size(); ++q)
200 for (
unsigned int d = 0;
d < dim; ++
d)
201 for (
unsigned int e = 0;
e < spacedim; ++
e)
202 if (
std::fabs(jac_0[
d][
e] - jac[
d][
e]) > zero_tolerance_double)
203 jacobian_constant =
false;
204 if (!jacobian_constant)
210 bool cell_cartesian = jacobian_constant;
211 for (
unsigned int d = 0;
d < dim; ++
d)
212 for (
unsigned int e = 0;
e < dim; ++
e)
214 if (
std::fabs(jac_0[
d][
e]) > zero_tolerance_double)
216 cell_cartesian =
false;
223 else if (jacobian_constant)
246 template <
int dim,
int spacedim = dim,
typename Number =
double>
254 Number>::vectorized_value_type;
339 template <
typename ContainerType>
342 const ContainerType &cell_iterator_range,
343 const std::vector<std::vector<
Point<dim>>> &unit_points_vector,
352 template <
typename ContainerType>
355 const ContainerType &cell_iterator_range,
363 template <
typename Iterator>
374 template <
typename Iterator>
428 get_JxW(
const unsigned int offset)
const;
458 boost::signals2::connection
464 template <
bool is_face>
467 const unsigned int face_number)
const;
486 const unsigned int geometry_index)
const;
522 template <
typename NumberType>
555 const unsigned int n_q_points,
564 const unsigned int n_q_points,
573 const unsigned int n_q_points,
576 const std::vector<double> &weights,
577 const unsigned int compressed_unit_point_index_offset,
578 const bool affine_cell);
589 template <
typename ContainerType,
typename QuadratureType>
592 const ContainerType &cell_iterator_range,
593 const std::vector<QuadratureType> &quadrature_vector,
594 const unsigned int n_unfiltered_cells,
597 const QuadratureType &quadrature,
598 MappingData &mapping_data)> &compute_mapping_data);
639 std::shared_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
668 std::vector<::internal::MatrixFreeFunctions::GeometryType>
cell_type;
764 template <
int dim,
int spacedim,
typename Number>
770 , update_flags(update_flags)
772 , additional_data(additional_data)
798 std::make_unique<typename MappingQ<dim, spacedim>::InternalData>(
799 mapping_q->get_degree());
805 template <
int dim,
int spacedim,
typename Number>
809 n_q_points_unvectorized.clear();
810 unit_points_index.clear();
811 data_index_offsets.clear();
812 compressed_data_index_offsets.clear();
818 template <
int dim,
int spacedim,
typename Number>
822 const std::vector<
Point<dim>> &unit_points_in)
829 template <
int dim,
int spacedim,
typename Number>
836 std::vector<
Point<dim>>(unit_points_in.begin(),
837 unit_points_in.end()));
842 template <
int dim,
int spacedim,
typename Number>
848 n_q_points_unvectorized.resize(1);
849 n_q_points_unvectorized[0] = quadrature.
size();
851 const unsigned int n_q_points =
852 compute_n_q_points<VectorizedArrayType>(n_q_points_unvectorized[0]);
854 const unsigned int n_q_points_data =
855 compute_n_q_points<Number>(n_q_points_unvectorized[0]);
858 resize_unit_points(n_q_points);
859 resize_data_fields(n_q_points_data);
864 n_q_points_unvectorized[0],
872 update_flags_mapping,
876 internal_mapping_data,
880 if (!quadrature.
empty() &&
895 n_q_points_unvectorized[0],
902 unit_points_index.push_back(0);
903 data_index_offsets.push_back(0);
904 compressed_data_index_offsets.push_back(0);
906 state = State::single_cell;
912 template <
int dim,
int spacedim,
typename Number>
913 template <
typename ContainerType>
916 const ContainerType &cell_iterator_range,
917 const std::vector<std::vector<
Point<dim>>> &unit_points_vector,
918 const unsigned int n_unfiltered_cells)
920 const unsigned int n_cells = unit_points_vector.size();
922 std::distance(cell_iterator_range.begin(),
923 cell_iterator_range.end()));
925 std::vector<Quadrature<dim>> quadrature_vector(
n_cells);
930 reinit_cells(cell_iterator_range, quadrature_vector, n_unfiltered_cells);
935 template <
int dim,
int spacedim,
typename Number>
936 template <
typename ContainerType,
typename QuadratureType>
939 const ContainerType &cell_iterator_range,
940 const std::vector<QuadratureType> &quadrature_vector,
941 const unsigned int n_unfiltered_cells,
944 const QuadratureType &quadrature,
949 do_cell_index_compression =
952 const unsigned int n_cells = quadrature_vector.size();
954 std::distance(cell_iterator_range.begin(),
955 cell_iterator_range.end()));
957 n_q_points_unvectorized.reserve(
n_cells);
961 if (additional_data.store_cells)
962 cell_level_and_indices.resize(
n_cells);
965 unit_points_index.reserve(
n_cells + 1);
966 unit_points_index.push_back(0);
967 data_index_offsets.reserve(
n_cells + 1);
968 data_index_offsets.push_back(0);
969 for (
const auto &quadrature : quadrature_vector)
971 const unsigned int n_points = quadrature.size();
972 n_q_points_unvectorized.push_back(n_points);
974 const unsigned int n_q_points =
975 compute_n_q_points<VectorizedArrayType>(n_points);
976 unit_points_index.push_back(unit_points_index.back() + n_q_points);
978 const unsigned int n_q_points_data =
979 compute_n_q_points<Number>(n_points);
980 data_index_offsets.push_back(data_index_offsets.back() +
984 const unsigned int n_unit_points = unit_points_index.back();
985 const unsigned int n_data_points = data_index_offsets.back();
988 resize_unit_points(n_unit_points);
989 resize_data_fields(n_data_points);
991 if (do_cell_index_compression)
992 cell_index_to_compressed_cell_index.resize(n_unfiltered_cells,
997 unsigned int size_compressed_data = 0;
999 for (
const auto &cell : cell_iterator_range)
1001 if (additional_data.store_cells)
1004 cell_level_and_indices[
cell_index] = {cell->level(), cell->index()};
1007 const auto &quadrature = quadrature_vector[
cell_index];
1008 const bool empty = quadrature.empty();
1011 const unsigned int n_q_points = compute_n_q_points<VectorizedArrayType>(
1013 store_unit_points(unit_points_index[
cell_index],
1016 quadrature.get_points());
1019 compute_mapping_data(cell, quadrature, mapping_data);
1022 const unsigned int n_q_points_data =
1023 compute_n_q_points<Number>(n_q_points_unvectorized[
cell_index]);
1029 cell_type.push_back(
1034 cell_type.push_back(
1040 const bool affine_cells =
1054 const auto comparison_result =
1064 comparison_result ==
1067 compressed_data_index_offsets.push_back(
1068 compressed_data_index_offsets.back());
1072 const unsigned int n_compressed_data_last_cell =
1076 compute_n_q_points<Number>(
1079 compressed_data_index_offsets.push_back(
1080 compressed_data_index_offsets.back() +
1081 n_compressed_data_last_cell);
1085 compressed_data_index_offsets.push_back(0);
1088 mapping_data_last_cell = mapping_data;
1090 store_mapping_data(data_index_offsets[
cell_index],
1094 quadrature.get_weights(),
1103 size_compressed_data = compressed_data_index_offsets.back() + 1;
1105 size_compressed_data =
1107 compressed_data_index_offsets.back() + n_q_points_data);
1109 if (do_cell_index_compression)
1110 cell_index_to_compressed_cell_index[cell->active_cell_index()] =
1118 jacobians.resize(size_compressed_data);
1119 jacobians.shrink_to_fit();
1123 inverse_jacobians.resize(size_compressed_data);
1124 inverse_jacobians.shrink_to_fit();
1127 state = State::cell_vector;
1133 template <
int dim,
int spacedim,
typename Number>
1134 template <
typename ContainerType>
1137 const ContainerType &cell_iterator_range,
1139 const unsigned int n_unfiltered_cells)
1141 auto compute_mapping_data_for_cells =
1149 update_flags_mapping,
1153 internal_mapping_data,
1157 do_reinit_cells<ContainerType, Quadrature<dim>>(
1158 cell_iterator_range,
1161 compute_mapping_data_for_cells);
1166 template <
int dim,
int spacedim,
typename Number>
1167 template <
typename Iterator>
1172 const unsigned int n_unfiltered_cells)
1175 additional_data.use_global_weights ==
false,
1177 "There is no known use-case for AdditionalData::use_global_weights=true and reinit_surface()"));
1184 auto compute_mapping_data_for_surface =
1191 update_flags_mapping,
1194 internal_mapping_data,
1199 cell_iterator_range,
1202 compute_mapping_data_for_surface);
1207 template <
int dim,
int spacedim,
typename Number>
1208 template <
typename Iterator>
1213 const unsigned int n_unfiltered_cells)
1219 do_cell_index_compression =
1222 const unsigned int n_cells = quadrature_vector.size();
1224 std::distance(cell_iterator_range.
begin(),
1225 cell_iterator_range.
end()));
1228 cell_index_offset.resize(
n_cells);
1229 unsigned int n_faces = 0;
1231 for (
const auto &cell : cell_iterator_range)
1234 n_faces += cell->n_faces();
1238 n_q_points_unvectorized.reserve(n_faces);
1240 cell_type.reserve(n_faces);
1243 unit_points_index.resize(n_faces + 1);
1244 data_index_offsets.resize(n_faces + 1);
1246 unsigned int n_unit_points = 0;
1247 unsigned int n_data_points = 0;
1248 for (
const auto &cell : cell_iterator_range)
1250 for (
const auto &f : cell->face_indices())
1252 const unsigned int current_face_index =
1255 unit_points_index[current_face_index] = n_unit_points;
1256 data_index_offsets[current_face_index] = n_data_points;
1258 const unsigned int n_points =
1260 n_q_points_unvectorized.push_back(n_points);
1262 const unsigned int n_q_points =
1263 compute_n_q_points<VectorizedArrayType>(n_points);
1264 n_unit_points += n_q_points;
1266 const unsigned int n_q_points_data =
1267 compute_n_q_points<Number>(n_points);
1268 n_data_points += n_q_points_data;
1273 unit_points_index[n_faces] = n_unit_points;
1274 data_index_offsets[n_faces] = n_data_points;
1277 if (do_cell_index_compression)
1278 cell_index_to_compressed_cell_index.resize(n_unfiltered_cells,
1283 resize_unit_points_faces(n_unit_points);
1284 resize_unit_points(n_unit_points);
1285 resize_data_fields(n_data_points);
1290 for (
const auto &cell : cell_iterator_range)
1292 const auto &quadratures_on_faces = quadrature_vector[
cell_index];
1294 Assert(quadratures_on_faces.size() == cell->n_faces(),
1298 for (
const auto &f : cell->face_indices())
1300 const auto &quadrature_on_face = quadratures_on_faces[f];
1302 const auto quadrature_on_cell =
1307 const unsigned int current_face_index =
1311 const unsigned int n_q_points =
1312 compute_n_q_points<VectorizedArrayType>(
1313 n_q_points_unvectorized[current_face_index]);
1314 store_unit_points(unit_points_index[current_face_index],
1316 n_q_points_unvectorized[current_face_index],
1317 quadrature_on_cell.get_points());
1319 store_unit_points_faces(unit_points_index[current_face_index],
1321 n_q_points_unvectorized[current_face_index],
1322 quadrature_on_face.get_points());
1326 update_flags_mapping,
1330 internal_mapping_data,
1333 cell_type.push_back(
1336 compressed_data_index_offsets.push_back(
1337 data_index_offsets[current_face_index]);
1339 const unsigned int n_q_points_data = compute_n_q_points<Number>(
1340 n_q_points_unvectorized[current_face_index]);
1341 store_mapping_data(data_index_offsets[current_face_index],
1343 n_q_points_unvectorized[current_face_index],
1345 quadrature_on_face.get_weights(),
1346 data_index_offsets[current_face_index],
1349 if (do_cell_index_compression)
1350 cell_index_to_compressed_cell_index[cell->active_cell_index()] =
1356 state = State::faces_on_cells_in_vector;
1362 template <
int dim,
int spacedim,
typename Number>
1366 return state == State::faces_on_cells_in_vector;
1371 template <
int dim,
int spacedim,
typename Number>
1374 const unsigned int geometry_index)
const
1376 return n_q_points_unvectorized[geometry_index];
1380 template <
int dim,
int spacedim,
typename Number>
1383 const unsigned int geometry_index)
const
1385 return cell_type[geometry_index];
1390 template <
int dim,
int spacedim,
typename Number>
1396 additional_data.store_cells,
1398 "Cells have been not stored. You can enable this by Additional::store_cells."));
1406 template <
int dim,
int spacedim,
typename Number>
1407 template <
typename NumberType>
1410 const unsigned int n_q_points_unvectorized)
1412 const unsigned int n_lanes =
1414 const unsigned int n_filled_lanes_last_batch =
1415 n_q_points_unvectorized % n_lanes;
1416 unsigned int n_q_points = n_q_points_unvectorized / n_lanes;
1417 if (n_filled_lanes_last_batch > 0)
1424 template <
int dim,
int spacedim,
typename Number>
1425 template <
bool is_face>
1429 const unsigned int face_number)
const
1434 const unsigned int compressed_cell_index =
1438 Assert(state == State::cell_vector,
1440 "This mapping info is not reinitialized for a cell vector!"));
1441 return compressed_cell_index;
1447 "cell_index has to be set if face_number is specified!"));
1448 Assert(state == State::faces_on_cells_in_vector,
1449 ExcMessage(
"This mapping info is not reinitialized for faces"
1450 " on cells in a vector!"));
1451 return cell_index_offset[compressed_cell_index] + face_number;
1457 template <
int dim,
int spacedim,
typename Number>
1462 if (do_cell_index_compression)
1466 ExcMessage(
"Mapping info object was not initialized for this"
1467 " active cell index!"));
1468 return cell_index_to_compressed_cell_index[
cell_index];
1475 template <
int dim,
int spacedim,
typename Number>
1478 const unsigned int unit_points_index_offset,
1479 const unsigned int n_q_points,
1480 const unsigned int n_q_points_unvectorized,
1483 const unsigned int n_lanes =
1486 for (
unsigned int q = 0; q < n_q_points; ++q)
1488 const unsigned int offset = unit_points_index_offset + q;
1489 for (
unsigned int v = 0;
1490 v < n_lanes && q * n_lanes + v < n_q_points_unvectorized;
1492 for (
unsigned int d = 0;
d < dim; ++
d)
1494 unit_points[offset][
d], v) = points[q * n_lanes + v][
d];
1500 template <
int dim,
int spacedim,
typename Number>
1503 const unsigned int unit_points_index_offset,
1504 const unsigned int n_q_points,
1505 const unsigned int n_q_points_unvectorized,
1508 const unsigned int n_lanes =
1511 for (
unsigned int q = 0; q < n_q_points; ++q)
1513 const unsigned int offset = unit_points_index_offset + q;
1514 for (
unsigned int v = 0;
1515 v < n_lanes && q * n_lanes + v < n_q_points_unvectorized;
1517 for (
unsigned int d = 0;
d < dim - 1; ++
d)
1519 unit_points_faces[offset][
d], v) = points[q * n_lanes + v][
d];
1525 template <
int dim,
int spacedim,
typename Number>
1528 const unsigned int unit_points_index_offset,
1529 const unsigned int n_q_points,
1530 const unsigned int n_q_points_unvectorized,
1532 const std::vector<double> &weights,
1533 const unsigned int compressed_unit_point_index_offset,
1534 const bool affine_cell)
1536 const unsigned int n_lanes =
1539 for (
unsigned int q = 0; q < n_q_points; ++q)
1541 const unsigned int offset = unit_points_index_offset + q;
1542 const unsigned int compressed_offset =
1543 compressed_unit_point_index_offset + q;
1544 for (
unsigned int v = 0;
1545 v < n_lanes && q * n_lanes + v < n_q_points_unvectorized;
1548 if (q == 0 || !affine_cell)
1551 for (
unsigned int d = 0;
d < dim; ++
d)
1552 for (
unsigned int s = 0; s < spacedim; ++s)
1554 jacobians[compressed_offset][
d][s], v) =
1555 mapping_data.
jacobians[q * n_lanes + v][
d][s];
1556 if (update_flags_mapping &
1558 for (
unsigned int d = 0;
d < dim; ++
d)
1559 for (
unsigned int s = 0; s < spacedim; ++s)
1561 inverse_jacobians[compressed_offset][s][
d], v) =
1566 if (additional_data.use_global_weights)
1569 JxW_values[offset], v) = weights[q * n_lanes + v];
1574 JxW_values[offset], v) =
1579 for (
unsigned int s = 0; s < spacedim; ++s)
1581 normal_vectors[offset][s], v) =
1584 for (
unsigned int s = 0; s < spacedim; ++s)
1586 real_points[offset][s], v) =
1594 template <
int dim,
int spacedim,
typename Number>
1597 const unsigned int n_unit_point_batches)
1599 unit_points.resize(n_unit_point_batches);
1604 template <
int dim,
int spacedim,
typename Number>
1607 const unsigned int n_unit_point_batches)
1609 unit_points_faces.resize(n_unit_point_batches);
1614 template <
int dim,
int spacedim,
typename Number>
1617 const unsigned int n_data_point_batches)
1620 jacobians.resize(n_data_point_batches);
1622 inverse_jacobians.resize(n_data_point_batches);
1624 JxW_values.resize(n_data_point_batches);
1626 normal_vectors.resize(n_data_point_batches);
1628 real_points.resize(n_data_point_batches);
1633 template <
int dim,
int spacedim,
typename Number>
1638 const unsigned int offset)
const
1640 return unit_points.data() + offset;
1645 template <
int dim,
int spacedim,
typename Number>
1650 const unsigned int offset)
const
1652 return unit_points_faces.data() + offset;
1657 template <
int dim,
int spacedim,
typename Number>
1660 const unsigned int offset)
const
1662 return real_points.data() + offset;
1667 template <
int dim,
int spacedim,
typename Number>
1670 const unsigned int geometry_index)
const
1672 return unit_points_index[geometry_index];
1677 template <
int dim,
int spacedim,
typename Number>
1680 const unsigned int geometry_index)
const
1682 return data_index_offsets[geometry_index];
1686 template <
int dim,
int spacedim,
typename Number>
1689 const unsigned int geometry_index)
const
1691 return compressed_data_index_offsets[geometry_index];
1696 template <
int dim,
int spacedim,
typename Number>
1699 const unsigned int offset)
const
1701 return jacobians.data() + offset;
1706 template <
int dim,
int spacedim,
typename Number>
1709 const unsigned int offset)
const
1711 return inverse_jacobians.data() + offset;
1715 template <
int dim,
int spacedim,
typename Number>
1718 const unsigned int offset)
const
1720 return normal_vectors.data() + offset;
1725 template <
int dim,
int spacedim,
typename Number>
1726 inline const Number *
1729 return JxW_values.data() + offset;
1734 template <
int dim,
int spacedim,
typename Number>
1743 template <
int dim,
int spacedim,
typename Number>
1747 return update_flags;
1752 template <
int dim,
int spacedim,
typename Number>
1756 return update_flags_mapping;
1761 template <
int dim,
int spacedim,
typename Number>
1762 boost::signals2::connection
1764 const std::function<
void()> &set_is_reinitialized)
1766 return is_reinitialized.connect(set_is_reinitialized);
1771 template <
int dim,
int spacedim,
typename Number>
1778 memory += cell_type.capacity() *
1791 cell_index_to_compressed_cell_index);
1793 memory +=
sizeof(*this);
IteratorOverIterators end() const
IteratorOverIterators begin()
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
const UpdateFlags update_flags
unsigned int compute_compressed_data_index_offset(const unsigned int geometry_index) const
AlignedVector< Point< dim, VectorizedArrayType > > unit_points
const DerivativeForm< 1, spacedim, dim, Number > * get_inverse_jacobian(const unsigned int offset) const
AlignedVector< Number > JxW_values
::internal::MatrixFreeFunctions::GeometryType get_cell_type(const unsigned int geometry_index) const
void reinit_surface(const IteratorRange< Iterator > &cell_iterator_range, const std::vector< ImmersedSurfaceQuadrature< dim >> &quadrature_vector, const unsigned int n_unfiltered_cells=numbers::invalid_unsigned_int)
std::shared_ptr< typename Mapping< dim, spacedim >::InternalDataBase > internal_mapping_data
const Point< dim, Number > * get_real_point(const unsigned int offset) const
boost::signals2::connection connect_is_reinitialized(const std::function< void()> &set_is_reinitialized)
const Mapping< dim, spacedim > & get_mapping() const
AlignedVector< DerivativeForm< 1, dim, spacedim, Number > > jacobians
SmartPointer< const Triangulation< dim, spacedim > > triangulation
const Point< dim, VectorizedArrayType > * get_unit_point(const unsigned int offset) const
unsigned int get_n_q_points_unvectorized(const unsigned int geometry_index) const
unsigned int compute_compressed_cell_index(const unsigned int cell_index) const
const Point< dim - 1, VectorizedArrayType > * get_unit_point_faces(const unsigned int offset) const
std::vector< unsigned int > cell_index_offset
MappingInfo(const Mapping< dim > &mapping, const UpdateFlags update_flags, const AdditionalData additional_data=AdditionalData())
const SmartPointer< const Mapping< dim, spacedim > > mapping
AlignedVector< Tensor< 1, spacedim, Number > > normal_vectors
unsigned int compute_geometry_index_offset(const unsigned int cell_index, const unsigned int face_number) const
std::vector<::internal::MatrixFreeFunctions::GeometryType > cell_type
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const std::vector< Point< dim >> &unit_points)
const DerivativeForm< 1, dim, spacedim, Number > * get_jacobian(const unsigned int offset) const
boost::signals2::signal< void()> is_reinitialized
void reinit_faces(const IteratorRange< Iterator > &cell_iterator_range, const std::vector< std::vector< Quadrature< dim - 1 >>> &quadrature_vector, const unsigned int n_unfiltered_cells=numbers::invalid_unsigned_int)
void resize_unit_points_faces(const unsigned int n_unit_point_batches)
const Tensor< 1, spacedim, Number > * get_normal_vector(const unsigned int offset) const
typename ::internal::VectorizedArrayTrait< Number >::vectorized_value_type VectorizedArrayType
void resize_unit_points(const unsigned int n_unit_point_batches)
std::size_t memory_consumption() const
UpdateFlags get_update_flags() const
void resize_data_fields(const unsigned int n_data_point_batches)
unsigned int compute_unit_point_index_offset(const unsigned int geometry_index) const
std::vector< unsigned int > cell_index_to_compressed_cell_index
AlignedVector< Point< spacedim, Number > > real_points
void store_mapping_data(const unsigned int unit_points_index_offset, const unsigned int n_q_points, const unsigned int n_q_points_unvectorized, const MappingData &mapping_data, const std::vector< double > &weights, const unsigned int compressed_unit_point_index_offset, const bool affine_cell)
unsigned int compute_data_index_offset(const unsigned int geometry_index) const
AlignedVector< DerivativeForm< 1, spacedim, dim, Number > > inverse_jacobians
@ faces_on_cells_in_vector
void reinit_cells(const ContainerType &cell_iterator_range, const std::vector< std::vector< Point< dim >>> &unit_points_vector, const unsigned int n_unfiltered_cells=numbers::invalid_unsigned_int)
std::vector< unsigned int > unit_points_index
std::vector< std::pair< int, int > > cell_level_and_indices
std::vector< unsigned int > compressed_data_index_offsets
bool do_cell_index_compression
std::vector< unsigned int > data_index_offsets
bool is_face_state() const
std::vector< unsigned int > n_q_points_unvectorized
const Number * get_JxW(const unsigned int offset) const
const AdditionalData additional_data
AlignedVector< Point< dim - 1, VectorizedArrayType > > unit_points_faces
unsigned int compute_n_q_points(const unsigned int n_q_points_unvectorized)
void store_unit_points_faces(const unsigned int unit_points_index_offset, const unsigned int n_q_points, const unsigned int n_q_points_unvectorized, const std::vector< Point< dim - 1 >> &points)
void store_unit_points(const unsigned int unit_points_index_offset, const unsigned int n_q_points, const unsigned int n_q_points_unvectorized, const std::vector< Point< dim >> &points)
UpdateFlags update_flags_mapping
void do_reinit_cells(const ContainerType &cell_iterator_range, const std::vector< QuadratureType > &quadrature_vector, const unsigned int n_unfiltered_cells, const std::function< void(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const QuadratureType &quadrature, MappingData &mapping_data)> &compute_mapping_data)
UpdateFlags get_update_flags_mapping() const
Triangulation< dim, spacedim >::cell_iterator get_cell_iterator(const unsigned int cell_index) const
static UpdateFlags required_update_flags(const SmartPointer< const Mapping< dim, spacedim >> mapping, const UpdateFlags &update_flags)
static void compute_mapping_data_for_face_quadrature(const SmartPointer< const Mapping< dim, spacedim >> mapping, const UpdateFlags &update_flags_mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim - 1 > &quadrature, std::shared_ptr< typename Mapping< dim, spacedim >::InternalDataBase > internal_mapping_data, MappingData &mapping_data)
static void compute_mapping_data_for_quadrature(const SmartPointer< const Mapping< dim, spacedim >> mapping, const UpdateFlags &update_flags_mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, CellSimilarity::Similarity &cell_similarity, const Quadrature< dim > &quadrature, std::shared_ptr< typename Mapping< dim, spacedim >::InternalDataBase > internal_mapping_data, MappingData &mapping_data)
static void compute_mapping_data_for_immersed_surface_quadrature(const SmartPointer< const Mapping< dim, spacedim >> mapping, const UpdateFlags &update_flags_mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ImmersedSurfaceQuadrature< dim > &quadrature, std::shared_ptr< typename Mapping< dim, spacedim >::InternalDataBase > internal_mapping_data, MappingData &mapping_data)
static void project_to_face(const ReferenceCell &reference_cell, const SubQuadrature &quadrature, const unsigned int face_no, std::vector< Point< dim >> &q_points)
const std::vector< Point< dim > > & get_points() const
const std::vector< double > & get_weights() const
unsigned int size() const
Triangulation< dim, spacedim > & get_triangulation()
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcMessage(std::string arg1)
@ update_normal_vectors
Normal vectors.
@ update_JxW_values
Transformed quadrature weights.
@ update_covariant_transformation
Covariant transformation.
@ update_jacobians
Volume element.
@ update_inverse_jacobians
Volume element.
@ update_gradients
Shape function gradients.
@ update_quadrature_points
Transformed quadrature points.
@ update_default
No update.
Expression fabs(const Expression &x)
@ general
No special properties.
std::enable_if_t< std::is_fundamental_v< T >, std::size_t > memory_consumption(const T &t)
::internal::MatrixFreeFunctions::GeometryType compute_geometry_type(const double diameter, const std::vector< DerivativeForm< 1, dim, spacedim, double >> &inverse_jacobians)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > epsilon(const Tensor< 2, dim, Number > &Grad_u)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
static const unsigned int invalid_unsigned_int
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
ComparisonResult compare(const std::vector< T > &v1, const std::vector< T > &v2) const
AdditionalData(const bool use_global_weights=false, const bool store_cells=false)
static constexpr std::size_t width()
static value_type & get(value_type &value, unsigned int c)