Reference documentation for deal.II version GIT 5c75851574 2022-09-25 16:30:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
fe_evaluation_data.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2020 - 2022 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #ifndef dealii_matrix_free_fe_evaluation_data_h
18 #define dealii_matrix_free_fe_evaluation_data_h
19 
20 
21 #include <deal.II/base/config.h>
22 
31 #include <deal.II/base/tensor.h>
33 
37 
38 
40 
41 
42 
43 namespace internal
44 {
46 
49  std::string,
50  << "You are requesting information from an FEEvaluation/FEFaceEvaluation "
51  << "object for which this kind of information has not been computed. What "
52  << "information these objects compute is determined by the update_* flags you "
53  << "pass to MatrixFree::reinit() via MatrixFree::AdditionalData. Here, "
54  << "the operation you are attempting requires the <" << arg1
55  << "> flag to be set, but it was apparently not specified "
56  << "upon initialization.");
57 } // namespace internal
58 
59 // forward declarations
60 template <int dim,
61  int fe_degree,
62  int n_q_points_1d = fe_degree + 1,
63  int n_components_ = 1,
64  typename Number = double,
65  typename VectorizedArrayType = VectorizedArray<Number>>
66 class FEEvaluation;
67 
68 template <int dim,
69  int n_components_,
70  typename Number,
71  bool = false,
72  typename = VectorizedArray<Number>>
73 class FEEvaluationBase;
74 
75 namespace internal
76 {
77  namespace MatrixFreeFunctions
78  {
79  template <int, typename>
80  class MappingDataOnTheFly;
81  }
82 } // namespace internal
83 
84 
113 template <int dim, typename Number, bool is_face>
115 {
117  using MappingInfoStorageType = internal::MatrixFreeFunctions::
118  MappingInfoStorage<(is_face ? dim - 1 : dim), dim, Number>;
120 
121 public:
122  static constexpr unsigned int dimension = dim;
123 
124  using NumberType = Number;
125 
126  using ScalarNumber =
128 
129  static constexpr unsigned int n_lanes = sizeof(Number) / sizeof(ScalarNumber);
130 
138  FEEvaluationData(const ShapeInfoType &shape_info,
139  const bool is_interior_face = true);
140 
144  FEEvaluationData(const FEEvaluationData &other) = default;
145 
151 
155  virtual ~FEEvaluationData() = default;
156 
162  void
164  const unsigned int n_components);
165 
170  void
173 
187 
192  Number
193  JxW(const unsigned int q_point) const;
194 
200  quadrature_point(const unsigned int q) const;
201 
214  inverse_jacobian(const unsigned int q_point) const;
215 
229  get_normal_vector(const unsigned int q_point) const;
230 
245  const Number *
247 
256  Number *
258 
269  const Number *
270  begin_values() const;
271 
282  Number *
284 
296  const Number *
298 
310  Number *
312 
325  const Number *
326  begin_hessians() const;
327 
340  Number *
342 
358  unsigned int
360 
368  get_cell_type() const;
369 
373  const ShapeInfoType &
374  get_shape_info() const;
375 
380  get_dof_info() const;
381 
387  const std::vector<unsigned int> &
389 
393  unsigned int
395 
399  unsigned int
401 
406  unsigned int
408 
413  unsigned int
415 
419  unsigned int
421 
434  std::uint8_t
435  get_face_no(const unsigned int v = 0) const;
436 
446  unsigned int
448 
457  std::uint8_t
458  get_face_orientation(const unsigned int v = 0) const;
459 
469 
478  bool
480 
485  const std::array<unsigned int, n_lanes> &
486  get_cell_ids() const
487  {
488 // implemented inline to avoid compilation problems on Windows
489 #ifdef DEBUG
491 #endif
492  return cell_ids;
493  }
494 
499  const std::array<unsigned int, n_lanes> &
500  get_face_ids() const
501  {
502 // implemented inline to avoid compilation problems on Windows
503 #ifdef DEBUG
505 #endif
506  return face_ids;
507  }
508 
513  unsigned int
515  {
516 // implemented inline to avoid compilation problems on Windows
517 #ifdef DEBUG
519 #endif
520 
521  return cell;
522  }
523 
528  const std::array<unsigned int, n_lanes> &
530  {
531 // implemented inline to avoid compilation problems on Windows
532 #ifdef DEBUG
534 #endif
535 
536  if (!is_face || dof_access_index ==
538  return cell_ids;
539  else
540  return face_ids;
541  }
542 
550 
564  Number
566 
573  void
574  set_cell_data(AlignedVector<Number> &array, const Number &value) const;
575 
580  template <typename T>
581  std::array<T, Number::size()>
583  const AlignedVector<std::array<T, Number::size()>> &array) const;
584 
589  template <typename T>
590  void
591  set_cell_data(AlignedVector<std::array<T, Number::size()>> &array,
592  const std::array<T, Number::size()> & value) const;
593 
602  Number
604 
613  void
614  set_face_data(AlignedVector<Number> &array, const Number &value) const;
615 
620  template <typename T>
621  std::array<T, Number::size()>
623  const AlignedVector<std::array<T, Number::size()>> &array) const;
624 
629  template <typename T>
630  void
631  set_face_data(AlignedVector<std::array<T, Number::size()>> &array,
632  const std::array<T, Number::size()> & value) const;
633 
641  {
643  const DoFInfo * dof_info;
645  unsigned int active_fe_index;
646  unsigned int active_quad_index;
648  };
649 
650 protected:
657  FEEvaluationData(const InitializationData &initialization_data,
658  const bool is_interior_face,
659  const unsigned int quad_no,
660  const unsigned int first_selected_component);
661 
667  const std::shared_ptr<
669  & mapping_data,
670  const unsigned int n_fe_components,
671  const unsigned int first_selected_component);
672 
680 
687 
695 
701  const unsigned int quad_no;
702 
707  const unsigned int n_fe_components;
708 
713  const unsigned int first_selected_component;
714 
718  const unsigned int active_fe_index;
719 
724  const unsigned int active_quad_index;
725 
732 
736  const unsigned int n_quadrature_points;
737 
743 
757 
762  const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>>
764 
769  const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>>
771 
778  const Number *J_value;
779 
784 
789 
794 
802 
813  Number *values_dofs;
814 
824  Number *values_quad;
825 
839 
851  Number *gradients_quad;
852 
865 
875  Number *hessians_quad;
876 
882 
889 
896 
903 
910 
917 
924 
931 
936  unsigned int cell;
937 
945 
951 
960  std::array<std::uint8_t, n_lanes> face_numbers;
961 
970  std::array<std::uint8_t, n_lanes> face_orientations;
971 
979  unsigned int subface_index;
980 
988 
993  std::array<unsigned int, n_lanes> cell_ids;
994 
999  std::array<unsigned int, n_lanes> face_ids;
1000 
1006  std::shared_ptr<
1009 
1015 
1016  // Make FEEvaluation and FEEvaluationBase objects friends for access to
1017  // protected member mapped_geometry.
1018  template <int, int, typename, bool, typename>
1019  friend class FEEvaluationBase;
1020 
1021  template <int, int, int, int, typename, typename>
1022  friend class FEEvaluation;
1023 };
1024 
1025 
1030 template <int dim,
1031  typename Number,
1032  bool is_face,
1033  typename VectorizedArrayType = VectorizedArray<Number>>
1036 
1037 
1038 /*----------------------- Inline functions ----------------------------------*/
1039 
1040 #ifndef DOXYGEN
1041 
1042 template <int dim, typename Number, bool is_face>
1044  const ShapeInfoType &shape_info,
1045  const bool is_interior_face)
1046  : FEEvaluationData(
1047  InitializationData{&shape_info, nullptr, nullptr, 0, 0, nullptr},
1049  0,
1050  0)
1051 {}
1052 
1053 
1054 template <int dim, typename Number, bool is_face>
1056  const InitializationData &initialization_data,
1057  const bool is_interior_face,
1058  const unsigned int quad_no,
1059  const unsigned int first_selected_component)
1060  : data(initialization_data.shape_info)
1061  , dof_info(initialization_data.dof_info)
1062  , mapping_data(initialization_data.mapping_data)
1063  , quad_no(quad_no)
1064  , n_fe_components(dof_info != nullptr ? dof_info->start_components.back() : 0)
1066  , active_fe_index(initialization_data.active_fe_index)
1067  , active_quad_index(initialization_data.active_quad_index)
1068  , descriptor(initialization_data.descriptor)
1069  , n_quadrature_points(descriptor == nullptr ?
1070  (is_face ? data->n_q_points_face : data->n_q_points) :
1072  , quadrature_points(nullptr)
1073  , jacobian(nullptr)
1074  , jacobian_gradients(nullptr)
1076  , J_value(nullptr)
1077  , normal_vectors(nullptr)
1078  , normal_x_jacobian(nullptr)
1080  descriptor != nullptr ? descriptor->quadrature_weights.begin() : nullptr)
1081 # ifdef DEBUG
1082  , is_reinitialized(false)
1083  , dof_values_initialized(false)
1084  , values_quad_initialized(false)
1086  , hessians_quad_initialized(false)
1087  , values_quad_submitted(false)
1088  , gradients_quad_submitted(false)
1089 # endif
1092  , dof_access_index(
1093  is_face ?
1094  (is_interior_face ?
1095  internal::MatrixFreeFunctions::DoFInfo::dof_access_face_interior :
1096  internal::MatrixFreeFunctions::DoFInfo::dof_access_face_exterior) :
1097  internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
1098  , subface_index(0)
1099  , cell_type(internal::MatrixFreeFunctions::general)
1100  , divergence_is_requested(false)
1101 {}
1102 
1103 
1104 
1105 template <int dim, typename Number, bool is_face>
1107  const std::shared_ptr<
1109  & mapped_geometry,
1110  const unsigned int n_fe_components,
1111  const unsigned int first_selected_component)
1112  : data(nullptr)
1113  , dof_info(nullptr)
1114  , mapping_data(nullptr)
1120  , descriptor(nullptr)
1122  mapped_geometry->get_data_storage().descriptor[0].n_q_points)
1123  , quadrature_points(nullptr)
1124  , jacobian(nullptr)
1125  , jacobian_gradients(nullptr)
1127  , J_value(nullptr)
1128  , normal_vectors(nullptr)
1129  , normal_x_jacobian(nullptr)
1130  , quadrature_weights(nullptr)
1131  , cell(0)
1132  , cell_type(internal::MatrixFreeFunctions::general)
1133  , interior_face(true)
1134  , dof_access_index(internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
1136  , is_reinitialized(false)
1137  , divergence_is_requested(false)
1138 {
1139  mapping_data = &mapped_geometry->get_data_storage();
1140  jacobian = mapped_geometry->get_data_storage().jacobians[0].begin();
1141  J_value = mapped_geometry->get_data_storage().JxW_values.begin();
1143  mapped_geometry->get_data_storage().jacobian_gradients[0].begin();
1144  jacobian_gradients_non_inverse = mapped_geometry->get_data_storage()
1145  .jacobian_gradients_non_inverse[0]
1146  .begin();
1148  mapped_geometry->get_data_storage().quadrature_points.begin();
1149 }
1150 
1151 
1152 
1153 template <int dim, typename Number, bool is_face>
1156 {
1163 
1164  data = other.data;
1165  dof_info = other.dof_info;
1166  mapping_data = other.mapping_data;
1167  descriptor = other.descriptor;
1168  jacobian = nullptr;
1169  J_value = nullptr;
1170  normal_vectors = nullptr;
1171  normal_x_jacobian = nullptr;
1172  jacobian_gradients = nullptr;
1174  quadrature_points = nullptr;
1176 
1177 # ifdef DEBUG
1178  is_reinitialized = false;
1179  dof_values_initialized = false;
1180  values_quad_initialized = false;
1182  hessians_quad_initialized = false;
1183  values_quad_submitted = false;
1184  gradients_quad_submitted = false;
1185 # endif
1186 
1188  interior_face = other.is_interior_face();
1190  is_face ?
1191  (is_interior_face() ?
1195  face_numbers[0] = 0;
1196  face_orientations[0] = 0;
1197  subface_index = 0;
1199  divergence_is_requested = false;
1200 
1201  return *this;
1202 }
1203 
1204 
1205 
1206 template <int dim, typename Number, bool is_face>
1207 inline void
1210  const unsigned int n_components)
1211 {
1213 
1214  const unsigned int tensor_dofs_per_component =
1215  Utilities::fixed_power<dim>(data->data.front().fe_degree + 1);
1216  const unsigned int dofs_per_component = data->dofs_per_component_on_cell;
1217 
1218  const unsigned int size_scratch_data =
1219  std::max(tensor_dofs_per_component + 1, dofs_per_component) * n_components *
1220  3 +
1221  2 * n_quadrature_points;
1222  const unsigned int size_data_arrays =
1224  (n_components * ((dim * (dim + 1)) / 2 + 2 * dim + 2) *
1226 
1227  const unsigned int allocated_size = size_scratch_data + size_data_arrays;
1228 # ifdef DEBUG
1230  scratch_data_array->resize(allocated_size,
1231  Number(numbers::signaling_nan<ScalarNumber>()));
1232 # else
1233  scratch_data_array->resize_fast(allocated_size);
1234 # endif
1235  scratch_data.reinit(scratch_data_array->begin() + size_data_arrays,
1236  size_scratch_data);
1237 
1238  // set the pointers to the correct position in the data array
1244  gradients_quad =
1250  hessians_quad =
1252  n_components * (dofs_per_component + (2 * dim + 2) * n_quadrature_points);
1253 }
1254 
1255 
1256 
1257 template <int dim, typename Number, bool is_face>
1258 inline void
1261 {
1262  Assert(is_face == true,
1263  ExcMessage("Faces can only be set if the is_face template parameter "
1264  "is true"));
1265  face_numbers[0] =
1267  subface_index = is_interior_face() == true ?
1269  face.subface_index;
1270 
1271  // First check if interior or exterior cell has non-standard orientation
1272  // (i.e. the third bit is one or not). Then set zero if this cell has
1273  // standard-orientation else copy the first three bits
1274  // (which is equivalent to modulo 8). See also the documentation of
1275  // internal::MatrixFreeFunctions::FaceToCellTopology::face_orientation.
1276  face_orientations[0] = (is_interior_face() == (face.face_orientation >= 8)) ?
1277  (face.face_orientation % 8) :
1278  0;
1279 
1280  if (is_interior_face())
1281  cell_ids = face.cells_interior;
1282  else
1283  cell_ids = face.cells_exterior;
1284 }
1285 
1286 
1287 
1288 template <int dim, typename Number, bool is_face>
1291  const unsigned int q_point) const
1292 {
1294  Assert(normal_vectors != nullptr,
1296  "update_normal_vectors"));
1298  return normal_vectors[0];
1299  else
1300  return normal_vectors[q_point];
1301 }
1302 
1303 
1304 
1305 template <int dim, typename Number, bool is_face>
1306 inline DEAL_II_ALWAYS_INLINE Number
1307 FEEvaluationData<dim, Number, is_face>::JxW(const unsigned int q_point) const
1308 {
1310  Assert(J_value != nullptr,
1312  "update_values|update_gradients"));
1314  {
1316  return J_value[0] * quadrature_weights[q_point];
1317  }
1318  else
1319  return J_value[q_point];
1320 }
1321 
1322 
1323 
1324 template <int dim, typename Number, bool is_face>
1327  const unsigned int q) const
1328 {
1330  Assert(this->quadrature_points != nullptr,
1332  "update_quadrature_points"));
1333 
1334  // Cartesian/affine mesh: only first vertex of cell is stored, we must
1335  // compute it through the Jacobian (which is stored in non-inverted and
1336  // non-transposed form as index '1' in the jacobian field)
1337  if (is_face == false &&
1339  {
1340  Assert(this->jacobian != nullptr, ExcNotInitialized());
1342 
1343  const Tensor<2, dim, Number> &jac = this->jacobian[1];
1345  for (unsigned int d = 0; d < dim; ++d)
1346  point[d] += jac[d][d] * static_cast<typename Number::value_type>(
1347  this->descriptor->quadrature.point(q)[d]);
1348  else
1349  for (unsigned int d = 0; d < dim; ++d)
1350  for (unsigned int e = 0; e < dim; ++e)
1351  point[d] += jac[d][e] * static_cast<typename Number::value_type>(
1352  this->descriptor->quadrature.point(q)[e]);
1353  return point;
1354  }
1355  else
1356  return this->quadrature_points[q];
1357 }
1358 
1359 
1360 
1361 template <int dim, typename Number, bool is_face>
1364  const unsigned int q_point) const
1365 {
1367  Assert(jacobian != nullptr,
1369  "update_gradients"));
1371  return jacobian[0];
1372  else
1373  return jacobian[q_point];
1374 }
1375 
1376 
1377 
1378 template <int dim, typename Number, bool is_face>
1379 inline const Number *
1381 {
1382  return values_dofs;
1383 }
1384 
1385 
1386 
1387 template <int dim, typename Number, bool is_face>
1388 inline Number *
1390 {
1391 # ifdef DEBUG
1392  dof_values_initialized = true;
1393 # endif
1394  return values_dofs;
1395 }
1396 
1397 
1398 
1399 template <int dim, typename Number, bool is_face>
1400 inline const Number *
1402 {
1403 # ifdef DEBUG
1405 # endif
1406  return values_quad;
1407 }
1408 
1409 
1410 
1411 template <int dim, typename Number, bool is_face>
1412 inline Number *
1414 {
1415 # ifdef DEBUG
1416  values_quad_initialized = true;
1417  values_quad_submitted = true;
1418 # endif
1419  return values_quad;
1420 }
1421 
1422 
1423 
1424 template <int dim, typename Number, bool is_face>
1425 inline const Number *
1427 {
1428 # ifdef DEBUG
1430  ExcNotInitialized());
1431 # endif
1432  return gradients_quad;
1433 }
1434 
1435 
1436 
1437 template <int dim, typename Number, bool is_face>
1438 inline Number *
1440 {
1441 # ifdef DEBUG
1442  gradients_quad_submitted = true;
1444 # endif
1445  return gradients_quad;
1446 }
1447 
1448 
1449 
1450 template <int dim, typename Number, bool is_face>
1451 inline const Number *
1453 {
1454 # ifdef DEBUG
1456 # endif
1457  return hessians_quad;
1458 }
1459 
1460 
1461 
1462 template <int dim, typename Number, bool is_face>
1463 inline Number *
1465 {
1466 # ifdef DEBUG
1468 # endif
1469  return hessians_quad;
1470 }
1471 
1472 
1473 
1474 template <int dim, typename Number, bool is_face>
1475 inline unsigned int
1477 {
1478  Assert(mapping_data != nullptr, ExcInternalError());
1479 
1480  if (dof_info == nullptr)
1481  return 0;
1482  else
1483  {
1486  }
1487 }
1488 
1489 
1490 
1491 template <int dim, typename Number, bool is_face>
1494 {
1495 # ifdef DEBUG
1497 # endif
1498  return cell_type;
1499 }
1500 
1501 
1502 
1503 template <int dim, typename Number, bool is_face>
1506 {
1507  Assert(data != nullptr, ExcInternalError());
1508  return *data;
1509 }
1510 
1511 
1512 
1513 template <int dim, typename Number, bool is_face>
1516 {
1517  Assert(dof_info != nullptr,
1518  ExcMessage(
1519  "FEEvaluation was not initialized with a MatrixFree object!"));
1520  return *dof_info;
1521 }
1522 
1523 
1524 
1525 template <int dim, typename Number, bool is_face>
1526 inline const std::vector<unsigned int> &
1528 {
1529  return data->lexicographic_numbering;
1530 }
1531 
1532 
1533 
1534 template <int dim, typename Number, bool is_face>
1535 inline unsigned int
1537 {
1538  return quad_no;
1539 }
1540 
1541 
1542 
1543 template <int dim, typename Number, bool is_face>
1544 inline unsigned int
1546 {
1547  if (is_face && dof_access_index ==
1550  else
1551  return cell;
1552 }
1553 
1554 
1555 
1556 template <int dim, typename Number, bool is_face>
1557 inline unsigned int
1559 {
1560  return active_fe_index;
1561 }
1562 
1563 
1564 
1565 template <int dim, typename Number, bool is_face>
1566 inline unsigned int
1568 {
1569  return active_quad_index;
1570 }
1571 
1572 
1573 
1574 template <int dim, typename Number, bool is_face>
1575 inline unsigned int
1577 {
1578  return first_selected_component;
1579 }
1580 
1581 
1582 
1583 template <int dim, typename Number, bool is_face>
1584 inline ArrayView<Number>
1586 {
1587  return scratch_data;
1588 }
1589 
1590 
1591 
1592 template <int dim, typename Number, bool is_face>
1593 std::uint8_t
1594 FEEvaluationData<dim, Number, is_face>::get_face_no(const unsigned int v) const
1595 {
1596  Assert(is_face, ExcNotInitialized());
1597  Assert(v == 0 || (cell != numbers::invalid_unsigned_int &&
1598  dof_access_index ==
1600  is_interior_face() == false),
1601  ExcMessage("All face numbers can only be queried for ECL at exterior "
1602  "faces. Use get_face_no() in other cases."));
1603 
1604  return face_numbers[v];
1605 }
1606 
1607 
1608 
1609 template <int dim, typename Number, bool is_face>
1610 inline unsigned int
1612 {
1613  return subface_index;
1614 }
1615 
1616 
1617 
1618 template <int dim, typename Number, bool is_face>
1619 std::uint8_t
1621  const unsigned int v) const
1622 {
1623  Assert(is_face, ExcNotInitialized());
1624  Assert(v == 0 || (cell != numbers::invalid_unsigned_int &&
1625  dof_access_index ==
1627  is_interior_face() == false),
1628  ExcMessage("All face numbers can only be queried for ECL at exterior "
1629  "faces. Use get_face_no() in other cases."));
1630 
1631  return face_orientations[v];
1632 }
1633 
1634 
1635 
1636 template <int dim, typename Number, bool is_face>
1639 {
1640  return dof_access_index;
1641 }
1642 
1643 
1644 
1645 template <int dim, typename Number, bool is_face>
1646 inline bool
1648 {
1649  return interior_face;
1650 }
1651 
1652 
1653 
1654 template <int dim, typename Number, bool is_face>
1657 {
1658  return {0U, n_quadrature_points};
1659 }
1660 
1661 
1662 
1663 namespace internal
1664 {
1665  template <std::size_t N,
1666  typename VectorizedArrayType2,
1667  typename GlobalVectorType,
1668  typename FU>
1669  inline void
1670  process_cell_or_face_data(const std::array<unsigned int, N> indices,
1671  GlobalVectorType & array,
1672  VectorizedArrayType2 & out,
1673  const FU & fu)
1674  {
1675  for (unsigned int i = 0; i < N; ++i)
1676  if (indices[i] != numbers::invalid_unsigned_int)
1677  {
1678  AssertIndexRange(indices[i] / N, array.size());
1679  fu(out[i], array[indices[i] / N][indices[i] % N]);
1680  }
1681  }
1682 } // namespace internal
1683 
1684 
1685 
1686 template <int dim, typename Number, bool is_face>
1687 inline Number
1689  const AlignedVector<Number> &array) const
1690 {
1691  Number out = Number(1.);
1692  internal::process_cell_or_face_data(this->get_cell_ids(),
1693  array,
1694  out,
1695  [](auto &local, const auto &global) {
1696  local = global;
1697  });
1698  return out;
1699 }
1700 
1701 
1702 
1703 template <int dim, typename Number, bool is_face>
1704 inline void
1706  AlignedVector<Number> &array,
1707  const Number & in) const
1708 {
1709  internal::process_cell_or_face_data(this->get_cell_ids(),
1710  array,
1711  in,
1712  [](const auto &local, auto &global) {
1713  global = local;
1714  });
1715 }
1716 
1717 
1718 
1719 template <int dim, typename Number, bool is_face>
1720 template <typename T>
1721 inline std::array<T, Number::size()>
1723  const AlignedVector<std::array<T, Number::size()>> &array) const
1724 {
1725  std::array<T, Number::size()> out;
1726  internal::process_cell_or_face_data(this->get_cell_ids(),
1727  array,
1728  out,
1729  [](auto &local, const auto &global) {
1730  local = global;
1731  });
1732  return out;
1733 }
1734 
1735 
1736 
1737 template <int dim, typename Number, bool is_face>
1738 template <typename T>
1739 inline void
1741  AlignedVector<std::array<T, Number::size()>> &array,
1742  const std::array<T, Number::size()> & in) const
1743 {
1744  internal::process_cell_or_face_data(this->get_cell_ids(),
1745  array,
1746  in,
1747  [](const auto &local, auto &global) {
1748  global = local;
1749  });
1750 }
1751 
1752 
1753 
1754 template <int dim, typename Number, bool is_face>
1755 inline Number
1757  const AlignedVector<Number> &array) const
1758 {
1759  Number out = Number(1.);
1760  internal::process_cell_or_face_data(this->get_face_ids(),
1761  array,
1762  out,
1763  [](auto &local, const auto &global) {
1764  local = global;
1765  });
1766  return out;
1767 }
1768 
1769 
1770 
1771 template <int dim, typename Number, bool is_face>
1772 inline void
1774  AlignedVector<Number> &array,
1775  const Number & in) const
1776 {
1777  internal::process_cell_or_face_data(this->get_face_ids(),
1778  array,
1779  in,
1780  [](const auto &local, auto &global) {
1781  global = local;
1782  });
1783 }
1784 
1785 
1786 
1787 template <int dim, typename Number, bool is_face>
1788 template <typename T>
1789 inline std::array<T, Number::size()>
1791  const AlignedVector<std::array<T, Number::size()>> &array) const
1792 {
1793  std::array<T, Number::size()> out;
1794  internal::process_cell_or_face_data(this->get_face_ids(),
1795  array,
1796  out,
1797  [](auto &local, const auto &global) {
1798  local = global;
1799  });
1800  return out;
1801 }
1802 
1803 
1804 
1805 template <int dim, typename Number, bool is_face>
1806 template <typename T>
1807 inline void
1809  AlignedVector<std::array<T, Number::size()>> &array,
1810  const std::array<T, Number::size()> & in) const
1811 {
1812  internal::process_cell_or_face_data(this->get_face_ids(),
1813  array,
1814  in,
1815  [](const auto &local, auto &global) {
1816  global = local;
1817  });
1818 }
1819 
1820 
1821 
1822 #endif // ifndef DOXYGEN
1823 
1824 
1826 
1827 #endif
void resize_fast(const size_type new_size)
iterator begin()
size_type size() const
void resize(const size_type new_size)
void reinit(value_type *starting_element, const std::size_t n_elements)
Definition: array_view.h:414
AlignedVector< VectorizedArrayType > * scratch_data_array
const unsigned int quad_no
const Number * J_value
Number * begin_values()
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
const Tensor< 1, dim *(dim+1)/2, Tensor< 1, dim, Number > > * jacobian_gradients_non_inverse
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
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
const DoFInfo * dof_info
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
Number * values_from_gradients_quad
const Number * begin_gradients() const
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
virtual ~FEEvaluationData()=default
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
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
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
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
Definition: config.h:102
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:457
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:458
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcAccessToUninitializedField()
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcNotInitialized()
#define Assert(cond, exc)
Definition: exceptions.h:1501
static ::ExceptionBase & ExcMatrixFreeAccessToUninitializedMappingField(std::string arg1)
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1695
#define AssertIndexRange(index, range)
Definition: exceptions.h:1760
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
static ::ExceptionBase & ExcMessage(std::string arg1)
static const char U
@ general
No special properties.
static const char T
static const char N
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:189
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={})
Definition: generators.cc:493
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
Definition: types.h:206
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
const MappingInfoStorageType * mapping_data
const MappingInfoStorageType::QuadratureDescriptor * descriptor
std::array< unsigned int, vectorization_width > cells_interior
Definition: face_info.h:61
std::array< unsigned int, vectorization_width > cells_exterior
Definition: face_info.h:76
std::vector< UnivariateShapeData< Number > > data
Definition: shape_info.h:429
std::vector< unsigned int > lexicographic_numbering
Definition: shape_info.h:423