Reference documentation for deal.II version GIT c415589cf0 2022-08-14 18:50:02+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 
157  void
159  const unsigned int n_components);
160 
165  void
168 
182 
187  Number
188  JxW(const unsigned int q_point) const;
189 
195  quadrature_point(const unsigned int q) const;
196 
209  inverse_jacobian(const unsigned int q_point) const;
210 
224  get_normal_vector(const unsigned int q_point) const;
225 
240  const Number *
242 
251  Number *
253 
264  const Number *
265  begin_values() const;
266 
277  Number *
279 
291  const Number *
293 
305  Number *
307 
320  const Number *
321  begin_hessians() const;
322 
335  Number *
337 
353  unsigned int
355 
363  get_cell_type() const;
364 
368  const ShapeInfoType &
369  get_shape_info() const;
370 
375  get_dof_info() const;
376 
382  const std::vector<unsigned int> &
384 
388  unsigned int
390 
394  unsigned int
396 
401  unsigned int
403 
408  unsigned int
410 
414  unsigned int
416 
429  std::uint8_t
430  get_face_no(const unsigned int v = 0) const;
431 
441  unsigned int
443 
452  std::uint8_t
453  get_face_orientation(const unsigned int v = 0) const;
454 
464 
473  bool
475 
480  const std::array<unsigned int, n_lanes> &
481  get_cell_ids() const
482  {
483 // implemented inline to avoid compilation problems on Windows
484 #ifdef DEBUG
486 #endif
487  return cell_ids;
488  }
489 
494  const std::array<unsigned int, n_lanes> &
495  get_face_ids() const
496  {
497 // implemented inline to avoid compilation problems on Windows
498 #ifdef DEBUG
500 #endif
501  return face_ids;
502  }
503 
508  unsigned int
510  {
511 // implemented inline to avoid compilation problems on Windows
512 #ifdef DEBUG
514 #endif
515 
516  return cell;
517  }
518 
523  const std::array<unsigned int, n_lanes> &
525  {
526 // implemented inline to avoid compilation problems on Windows
527 #ifdef DEBUG
529 #endif
530 
531  if (!is_face || dof_access_index ==
533  return cell_ids;
534  else
535  return face_ids;
536  }
537 
545 
559  Number
561 
568  void
569  set_cell_data(AlignedVector<Number> &array, const Number &value) const;
570 
575  template <typename T>
576  std::array<T, Number::size()>
578  const AlignedVector<std::array<T, Number::size()>> &array) const;
579 
584  template <typename T>
585  void
586  set_cell_data(AlignedVector<std::array<T, Number::size()>> &array,
587  const std::array<T, Number::size()> & value) const;
588 
597  Number
599 
608  void
609  set_face_data(AlignedVector<Number> &array, const Number &value) const;
610 
615  template <typename T>
616  std::array<T, Number::size()>
618  const AlignedVector<std::array<T, Number::size()>> &array) const;
619 
624  template <typename T>
625  void
626  set_face_data(AlignedVector<std::array<T, Number::size()>> &array,
627  const std::array<T, Number::size()> & value) const;
628 
636  {
638  const DoFInfo * dof_info;
640  unsigned int active_fe_index;
641  unsigned int active_quad_index;
643  };
644 
645 protected:
652  FEEvaluationData(const InitializationData &initialization_data,
653  const bool is_interior_face,
654  const unsigned int quad_no,
655  const unsigned int first_selected_component);
656 
662  const std::shared_ptr<
664  & mapping_data,
665  const unsigned int n_fe_components,
666  const unsigned int first_selected_component);
667 
675 
682 
690 
696  const unsigned int quad_no;
697 
702  const unsigned int n_fe_components;
703 
708  const unsigned int first_selected_component;
709 
713  const unsigned int active_fe_index;
714 
719  const unsigned int active_quad_index;
720 
727 
731  const unsigned int n_quadrature_points;
732 
738 
752 
757  const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>>
759 
764  const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>>
766 
773  const Number *J_value;
774 
779 
784 
789 
797 
808  Number *values_dofs;
809 
819  Number *values_quad;
820 
834 
846  Number *gradients_quad;
847 
860 
870  Number *hessians_quad;
871 
877 
884 
891 
898 
905 
912 
919 
926 
931  unsigned int cell;
932 
940 
946 
955  std::array<std::uint8_t, n_lanes> face_numbers;
956 
965  std::array<std::uint8_t, n_lanes> face_orientations;
966 
974  unsigned int subface_index;
975 
983 
988  std::array<unsigned int, n_lanes> cell_ids;
989 
994  std::array<unsigned int, n_lanes> face_ids;
995 
1001  std::shared_ptr<
1004 
1010 
1011  // Make FEEvaluation and FEEvaluationBase objects friends for access to
1012  // protected member mapped_geometry.
1013  template <int, int, typename, bool, typename>
1014  friend class FEEvaluationBase;
1015 
1016  template <int, int, int, int, typename, typename>
1017  friend class FEEvaluation;
1018 };
1019 
1020 
1025 template <int dim,
1026  typename Number,
1027  bool is_face,
1028  typename VectorizedArrayType = VectorizedArray<Number>>
1031 
1032 
1033 /*----------------------- Inline functions ----------------------------------*/
1034 
1035 #ifndef DOXYGEN
1036 
1037 template <int dim, typename Number, bool is_face>
1039  const ShapeInfoType &shape_info,
1040  const bool is_interior_face)
1041  : FEEvaluationData(
1042  InitializationData{&shape_info, nullptr, nullptr, 0, 0, nullptr},
1044  0,
1045  0)
1046 {}
1047 
1048 
1049 template <int dim, typename Number, bool is_face>
1051  const InitializationData &initialization_data,
1052  const bool is_interior_face,
1053  const unsigned int quad_no,
1054  const unsigned int first_selected_component)
1055  : data(initialization_data.shape_info)
1056  , dof_info(initialization_data.dof_info)
1057  , mapping_data(initialization_data.mapping_data)
1058  , quad_no(quad_no)
1059  , n_fe_components(dof_info != nullptr ? dof_info->start_components.back() : 0)
1061  , active_fe_index(initialization_data.active_fe_index)
1062  , active_quad_index(initialization_data.active_quad_index)
1063  , descriptor(initialization_data.descriptor)
1064  , n_quadrature_points(descriptor == nullptr ?
1065  (is_face ? data->n_q_points_face : data->n_q_points) :
1067  , quadrature_points(nullptr)
1068  , jacobian(nullptr)
1069  , jacobian_gradients(nullptr)
1071  , J_value(nullptr)
1072  , normal_vectors(nullptr)
1073  , normal_x_jacobian(nullptr)
1075  descriptor != nullptr ? descriptor->quadrature_weights.begin() : nullptr)
1076 # ifdef DEBUG
1077  , is_reinitialized(false)
1078  , dof_values_initialized(false)
1079  , values_quad_initialized(false)
1081  , hessians_quad_initialized(false)
1082  , values_quad_submitted(false)
1083  , gradients_quad_submitted(false)
1084 # endif
1087  , dof_access_index(
1088  is_face ?
1089  (is_interior_face ?
1090  internal::MatrixFreeFunctions::DoFInfo::dof_access_face_interior :
1091  internal::MatrixFreeFunctions::DoFInfo::dof_access_face_exterior) :
1092  internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
1093  , subface_index(0)
1094  , cell_type(internal::MatrixFreeFunctions::general)
1095  , divergence_is_requested(false)
1096 {}
1097 
1098 
1099 
1100 template <int dim, typename Number, bool is_face>
1102  const std::shared_ptr<
1104  & mapped_geometry,
1105  const unsigned int n_fe_components,
1106  const unsigned int first_selected_component)
1107  : data(nullptr)
1108  , dof_info(nullptr)
1109  , mapping_data(nullptr)
1115  , descriptor(nullptr)
1117  mapped_geometry->get_data_storage().descriptor[0].n_q_points)
1118  , quadrature_points(nullptr)
1119  , jacobian(nullptr)
1120  , jacobian_gradients(nullptr)
1122  , J_value(nullptr)
1123  , normal_vectors(nullptr)
1124  , normal_x_jacobian(nullptr)
1125  , quadrature_weights(nullptr)
1126  , cell(0)
1127  , cell_type(internal::MatrixFreeFunctions::general)
1128  , interior_face(true)
1129  , dof_access_index(internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
1131  , is_reinitialized(false)
1132  , divergence_is_requested(false)
1133 {
1134  mapping_data = &mapped_geometry->get_data_storage();
1135  jacobian = mapped_geometry->get_data_storage().jacobians[0].begin();
1136  J_value = mapped_geometry->get_data_storage().JxW_values.begin();
1138  mapped_geometry->get_data_storage().jacobian_gradients[0].begin();
1139  jacobian_gradients_non_inverse = mapped_geometry->get_data_storage()
1140  .jacobian_gradients_non_inverse[0]
1141  .begin();
1143  mapped_geometry->get_data_storage().quadrature_points.begin();
1144 }
1145 
1146 
1147 
1148 template <int dim, typename Number, bool is_face>
1151 {
1158 
1159  data = other.data;
1160  dof_info = other.dof_info;
1161  mapping_data = other.mapping_data;
1162  descriptor = other.descriptor;
1163  jacobian = nullptr;
1164  J_value = nullptr;
1165  normal_vectors = nullptr;
1166  normal_x_jacobian = nullptr;
1167  jacobian_gradients = nullptr;
1169  quadrature_points = nullptr;
1171 
1172 # ifdef DEBUG
1173  is_reinitialized = false;
1174  dof_values_initialized = false;
1175  values_quad_initialized = false;
1177  hessians_quad_initialized = false;
1178  values_quad_submitted = false;
1179  gradients_quad_submitted = false;
1180 # endif
1181 
1183  interior_face = other.is_interior_face();
1185  is_face ?
1186  (is_interior_face() ?
1190  face_numbers[0] = 0;
1191  face_orientations[0] = 0;
1192  subface_index = 0;
1194  divergence_is_requested = false;
1195 
1196  return *this;
1197 }
1198 
1199 
1200 
1201 template <int dim, typename Number, bool is_face>
1202 inline void
1205  const unsigned int n_components)
1206 {
1208 
1209  const unsigned int tensor_dofs_per_component =
1210  Utilities::fixed_power<dim>(data->data.front().fe_degree + 1);
1211  const unsigned int dofs_per_component = data->dofs_per_component_on_cell;
1212 
1213  const unsigned int size_scratch_data =
1214  std::max(tensor_dofs_per_component + 1, dofs_per_component) * n_components *
1215  3 +
1216  2 * n_quadrature_points;
1217  const unsigned int size_data_arrays =
1219  (n_components * ((dim * (dim + 1)) / 2 + 2 * dim + 2) *
1221 
1222  const unsigned int allocated_size = size_scratch_data + size_data_arrays;
1223 # ifdef DEBUG
1225  scratch_data_array->resize(allocated_size,
1226  Number(numbers::signaling_nan<ScalarNumber>()));
1227 # else
1228  scratch_data_array->resize_fast(allocated_size);
1229 # endif
1230  scratch_data.reinit(scratch_data_array->begin() + size_data_arrays,
1231  size_scratch_data);
1232 
1233  // set the pointers to the correct position in the data array
1239  gradients_quad =
1245  hessians_quad =
1247  n_components * (dofs_per_component + (2 * dim + 2) * n_quadrature_points);
1248 }
1249 
1250 
1251 
1252 template <int dim, typename Number, bool is_face>
1253 inline void
1256 {
1257  Assert(is_face == true,
1258  ExcMessage("Faces can only be set if the is_face template parameter "
1259  "is true"));
1260  face_numbers[0] =
1262  subface_index = is_interior_face() == true ?
1264  face.subface_index;
1265 
1266  // First check if interior or exterior cell has non-standard orientation
1267  // (i.e. the third bit is one or not). Then set zero if this cell has
1268  // standard-orientation else copy the first three bits
1269  // (which is equivalent to modulo 8). See also the documentation of
1270  // internal::MatrixFreeFunctions::FaceToCellTopology::face_orientation.
1271  face_orientations[0] = (is_interior_face() == (face.face_orientation >= 8)) ?
1272  (face.face_orientation % 8) :
1273  0;
1274 
1275  if (is_interior_face())
1276  cell_ids = face.cells_interior;
1277  else
1278  cell_ids = face.cells_exterior;
1279 }
1280 
1281 
1282 
1283 template <int dim, typename Number, bool is_face>
1286  const unsigned int q_point) const
1287 {
1289  Assert(normal_vectors != nullptr,
1291  "update_normal_vectors"));
1293  return normal_vectors[0];
1294  else
1295  return normal_vectors[q_point];
1296 }
1297 
1298 
1299 
1300 template <int dim, typename Number, bool is_face>
1301 inline DEAL_II_ALWAYS_INLINE Number
1302 FEEvaluationData<dim, Number, is_face>::JxW(const unsigned int q_point) const
1303 {
1305  Assert(J_value != nullptr,
1307  "update_values|update_gradients"));
1309  {
1311  return J_value[0] * quadrature_weights[q_point];
1312  }
1313  else
1314  return J_value[q_point];
1315 }
1316 
1317 
1318 
1319 template <int dim, typename Number, bool is_face>
1322  const unsigned int q) const
1323 {
1325  Assert(this->quadrature_points != nullptr,
1327  "update_quadrature_points"));
1328 
1329  // Cartesian/affine mesh: only first vertex of cell is stored, we must
1330  // compute it through the Jacobian (which is stored in non-inverted and
1331  // non-transposed form as index '1' in the jacobian field)
1332  if (is_face == false &&
1334  {
1335  Assert(this->jacobian != nullptr, ExcNotInitialized());
1337 
1338  const Tensor<2, dim, Number> &jac = this->jacobian[1];
1340  for (unsigned int d = 0; d < dim; ++d)
1341  point[d] += jac[d][d] * static_cast<typename Number::value_type>(
1342  this->descriptor->quadrature.point(q)[d]);
1343  else
1344  for (unsigned int d = 0; d < dim; ++d)
1345  for (unsigned int e = 0; e < dim; ++e)
1346  point[d] += jac[d][e] * static_cast<typename Number::value_type>(
1347  this->descriptor->quadrature.point(q)[e]);
1348  return point;
1349  }
1350  else
1351  return this->quadrature_points[q];
1352 }
1353 
1354 
1355 
1356 template <int dim, typename Number, bool is_face>
1359  const unsigned int q_point) const
1360 {
1362  Assert(jacobian != nullptr,
1364  "update_gradients"));
1366  return jacobian[0];
1367  else
1368  return jacobian[q_point];
1369 }
1370 
1371 
1372 
1373 template <int dim, typename Number, bool is_face>
1374 inline const Number *
1376 {
1377  return values_dofs;
1378 }
1379 
1380 
1381 
1382 template <int dim, typename Number, bool is_face>
1383 inline Number *
1385 {
1386 # ifdef DEBUG
1387  dof_values_initialized = true;
1388 # endif
1389  return values_dofs;
1390 }
1391 
1392 
1393 
1394 template <int dim, typename Number, bool is_face>
1395 inline const Number *
1397 {
1398 # ifdef DEBUG
1400 # endif
1401  return values_quad;
1402 }
1403 
1404 
1405 
1406 template <int dim, typename Number, bool is_face>
1407 inline Number *
1409 {
1410 # ifdef DEBUG
1411  values_quad_initialized = true;
1412  values_quad_submitted = true;
1413 # endif
1414  return values_quad;
1415 }
1416 
1417 
1418 
1419 template <int dim, typename Number, bool is_face>
1420 inline const Number *
1422 {
1423 # ifdef DEBUG
1425  ExcNotInitialized());
1426 # endif
1427  return gradients_quad;
1428 }
1429 
1430 
1431 
1432 template <int dim, typename Number, bool is_face>
1433 inline Number *
1435 {
1436 # ifdef DEBUG
1437  gradients_quad_submitted = true;
1439 # endif
1440  return gradients_quad;
1441 }
1442 
1443 
1444 
1445 template <int dim, typename Number, bool is_face>
1446 inline const Number *
1448 {
1449 # ifdef DEBUG
1451 # endif
1452  return hessians_quad;
1453 }
1454 
1455 
1456 
1457 template <int dim, typename Number, bool is_face>
1458 inline Number *
1460 {
1461 # ifdef DEBUG
1463 # endif
1464  return hessians_quad;
1465 }
1466 
1467 
1468 
1469 template <int dim, typename Number, bool is_face>
1470 inline unsigned int
1472 {
1473  Assert(mapping_data != nullptr, ExcInternalError());
1474 
1475  if (dof_info == nullptr)
1476  return 0;
1477  else
1478  {
1481  }
1482 }
1483 
1484 
1485 
1486 template <int dim, typename Number, bool is_face>
1489 {
1490 # ifdef DEBUG
1492 # endif
1493  return cell_type;
1494 }
1495 
1496 
1497 
1498 template <int dim, typename Number, bool is_face>
1501 {
1502  Assert(data != nullptr, ExcInternalError());
1503  return *data;
1504 }
1505 
1506 
1507 
1508 template <int dim, typename Number, bool is_face>
1511 {
1512  Assert(dof_info != nullptr,
1513  ExcMessage(
1514  "FEEvaluation was not initialized with a MatrixFree object!"));
1515  return *dof_info;
1516 }
1517 
1518 
1519 
1520 template <int dim, typename Number, bool is_face>
1521 inline const std::vector<unsigned int> &
1523 {
1524  return data->lexicographic_numbering;
1525 }
1526 
1527 
1528 
1529 template <int dim, typename Number, bool is_face>
1530 inline unsigned int
1532 {
1533  return quad_no;
1534 }
1535 
1536 
1537 
1538 template <int dim, typename Number, bool is_face>
1539 inline unsigned int
1541 {
1542  if (is_face && dof_access_index ==
1545  else
1546  return cell;
1547 }
1548 
1549 
1550 
1551 template <int dim, typename Number, bool is_face>
1552 inline unsigned int
1554 {
1555  return active_fe_index;
1556 }
1557 
1558 
1559 
1560 template <int dim, typename Number, bool is_face>
1561 inline unsigned int
1563 {
1564  return active_quad_index;
1565 }
1566 
1567 
1568 
1569 template <int dim, typename Number, bool is_face>
1570 inline unsigned int
1572 {
1573  return first_selected_component;
1574 }
1575 
1576 
1577 
1578 template <int dim, typename Number, bool is_face>
1579 inline ArrayView<Number>
1581 {
1582  return scratch_data;
1583 }
1584 
1585 
1586 
1587 template <int dim, typename Number, bool is_face>
1588 std::uint8_t
1589 FEEvaluationData<dim, Number, is_face>::get_face_no(const unsigned int v) const
1590 {
1591  Assert(is_face, ExcNotInitialized());
1592  Assert(v == 0 || (cell != numbers::invalid_unsigned_int &&
1593  dof_access_index ==
1595  is_interior_face() == false),
1596  ExcMessage("All face numbers can only be queried for ECL at exterior "
1597  "faces. Use get_face_no() in other cases."));
1598 
1599  return face_numbers[v];
1600 }
1601 
1602 
1603 
1604 template <int dim, typename Number, bool is_face>
1605 inline unsigned int
1607 {
1608  return subface_index;
1609 }
1610 
1611 
1612 
1613 template <int dim, typename Number, bool is_face>
1614 std::uint8_t
1616  const unsigned int v) const
1617 {
1618  Assert(is_face, ExcNotInitialized());
1619  Assert(v == 0 || (cell != numbers::invalid_unsigned_int &&
1620  dof_access_index ==
1622  is_interior_face() == false),
1623  ExcMessage("All face numbers can only be queried for ECL at exterior "
1624  "faces. Use get_face_no() in other cases."));
1625 
1626  return face_orientations[v];
1627 }
1628 
1629 
1630 
1631 template <int dim, typename Number, bool is_face>
1634 {
1635  return dof_access_index;
1636 }
1637 
1638 
1639 
1640 template <int dim, typename Number, bool is_face>
1641 inline bool
1643 {
1644  return interior_face;
1645 }
1646 
1647 
1648 
1649 template <int dim, typename Number, bool is_face>
1652 {
1653  return {0U, n_quadrature_points};
1654 }
1655 
1656 
1657 
1658 namespace internal
1659 {
1660  template <std::size_t N,
1661  typename VectorizedArrayType2,
1662  typename GlobalVectorType,
1663  typename FU>
1664  inline void
1665  process_cell_or_face_data(const std::array<unsigned int, N> indices,
1666  GlobalVectorType & array,
1667  VectorizedArrayType2 & out,
1668  const FU & fu)
1669  {
1670  for (unsigned int i = 0; i < N; ++i)
1671  if (indices[i] != numbers::invalid_unsigned_int)
1672  {
1673  AssertIndexRange(indices[i] / N, array.size());
1674  fu(out[i], array[indices[i] / N][indices[i] % N]);
1675  }
1676  }
1677 } // namespace internal
1678 
1679 
1680 
1681 template <int dim, typename Number, bool is_face>
1682 inline Number
1684  const AlignedVector<Number> &array) const
1685 {
1686  Number out = Number(1.);
1687  internal::process_cell_or_face_data(this->get_cell_ids(),
1688  array,
1689  out,
1690  [](auto &local, const auto &global) {
1691  local = global;
1692  });
1693  return out;
1694 }
1695 
1696 
1697 
1698 template <int dim, typename Number, bool is_face>
1699 inline void
1701  AlignedVector<Number> &array,
1702  const Number & in) const
1703 {
1704  internal::process_cell_or_face_data(this->get_cell_ids(),
1705  array,
1706  in,
1707  [](const auto &local, auto &global) {
1708  global = local;
1709  });
1710 }
1711 
1712 
1713 
1714 template <int dim, typename Number, bool is_face>
1715 template <typename T>
1716 inline std::array<T, Number::size()>
1718  const AlignedVector<std::array<T, Number::size()>> &array) const
1719 {
1720  std::array<T, Number::size()> out;
1721  internal::process_cell_or_face_data(this->get_cell_ids(),
1722  array,
1723  out,
1724  [](auto &local, const auto &global) {
1725  local = global;
1726  });
1727  return out;
1728 }
1729 
1730 
1731 
1732 template <int dim, typename Number, bool is_face>
1733 template <typename T>
1734 inline void
1736  AlignedVector<std::array<T, Number::size()>> &array,
1737  const std::array<T, Number::size()> & in) const
1738 {
1739  internal::process_cell_or_face_data(this->get_cell_ids(),
1740  array,
1741  in,
1742  [](const auto &local, auto &global) {
1743  global = local;
1744  });
1745 }
1746 
1747 
1748 
1749 template <int dim, typename Number, bool is_face>
1750 inline Number
1752  const AlignedVector<Number> &array) const
1753 {
1754  Number out = Number(1.);
1755  internal::process_cell_or_face_data(this->get_face_ids(),
1756  array,
1757  out,
1758  [](auto &local, const auto &global) {
1759  local = global;
1760  });
1761  return out;
1762 }
1763 
1764 
1765 
1766 template <int dim, typename Number, bool is_face>
1767 inline void
1769  AlignedVector<Number> &array,
1770  const Number & in) const
1771 {
1772  internal::process_cell_or_face_data(this->get_face_ids(),
1773  array,
1774  in,
1775  [](const auto &local, auto &global) {
1776  global = local;
1777  });
1778 }
1779 
1780 
1781 
1782 template <int dim, typename Number, bool is_face>
1783 template <typename T>
1784 inline std::array<T, Number::size()>
1786  const AlignedVector<std::array<T, Number::size()>> &array) const
1787 {
1788  std::array<T, Number::size()> out;
1789  internal::process_cell_or_face_data(this->get_face_ids(),
1790  array,
1791  out,
1792  [](auto &local, const auto &global) {
1793  local = global;
1794  });
1795  return out;
1796 }
1797 
1798 
1799 
1800 template <int dim, typename Number, bool is_face>
1801 template <typename T>
1802 inline void
1804  AlignedVector<std::array<T, Number::size()>> &array,
1805  const std::array<T, Number::size()> & in) const
1806 {
1807  internal::process_cell_or_face_data(this->get_face_ids(),
1808  array,
1809  in,
1810  [](const auto &local, auto &global) {
1811  global = local;
1812  });
1813 }
1814 
1815 
1816 
1817 #endif // ifndef DOXYGEN
1818 
1819 
1821 
1822 #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
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:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcAccessToUninitializedField()
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcNotInitialized()
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcMatrixFreeAccessToUninitializedMappingField(std::string arg1)
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1667
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
#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:190
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:201
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:423
std::vector< unsigned int > lexicographic_numbering
Definition: shape_info.h:417