Reference documentation for deal.II version GIT 1f418c7800 2022-06-28 06:15: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 
30 #include <deal.II/base/tensor.h>
32 
36 
37 
39 
40 
41 
42 namespace internal
43 {
45 
48  std::string,
49  << "You are requesting information from an FEEvaluation/FEFaceEvaluation "
50  << "object for which this kind of information has not been computed. What "
51  << "information these objects compute is determined by the update_* flags you "
52  << "pass to MatrixFree::reinit() via MatrixFree::AdditionalData. Here, "
53  << "the operation you are attempting requires the <" << arg1
54  << "> flag to be set, but it was apparently not specified "
55  << "upon initialization.");
56 } // namespace internal
57 
58 // forward declarations
59 template <int dim,
60  int fe_degree,
61  int n_q_points_1d = fe_degree + 1,
62  int n_components_ = 1,
63  typename Number = double,
64  typename VectorizedArrayType = VectorizedArray<Number>>
65 class FEEvaluation;
66 
67 template <int dim,
68  int n_components_,
69  typename Number,
70  bool = false,
71  typename = VectorizedArray<Number>>
72 class FEEvaluationBase;
73 
74 namespace internal
75 {
76  namespace MatrixFreeFunctions
77  {
78  template <int, typename>
79  class MappingDataOnTheFly;
80  }
81 } // namespace internal
82 
83 
112 template <int dim, typename Number, bool is_face>
114 {
116  using MappingInfoStorageType = internal::MatrixFreeFunctions::
117  MappingInfoStorage<(is_face ? dim - 1 : dim), dim, Number>;
119 
120 public:
121  static constexpr unsigned int dimension = dim;
122 
123  using NumberType = Number;
124 
125  using ScalarNumber =
127 
128  static constexpr unsigned int n_lanes = sizeof(Number) / sizeof(ScalarNumber);
129 
137  FEEvaluationData(const ShapeInfoType &shape_info,
138  const bool is_interior_face = true);
139 
143  FEEvaluationData(const FEEvaluationData &other) = default;
144 
150 
156  void
158  const unsigned int n_components);
159 
164  void
167 
172 
181 
186  Number
187  JxW(const unsigned int q_point) const;
188 
194  quadrature_point(const unsigned int q) const;
195 
208  inverse_jacobian(const unsigned int q_point) const;
209 
223  get_normal_vector(const unsigned int q_point) const;
224 
226 
239  const Number *
241 
250  Number *
252 
263  const Number *
264  begin_values() const;
265 
276  Number *
278 
290  const Number *
292 
304  Number *
306 
319  const Number *
320  begin_hessians() const;
321 
334  Number *
336 
338 
343 
352  unsigned int
354 
362  get_cell_type() const;
363 
367  const ShapeInfoType &
368  get_shape_info() const;
369 
374  get_dof_info() const;
375 
381  const std::vector<unsigned int> &
383 
387  unsigned int
389 
393  unsigned int
395 
400  unsigned int
402 
407  unsigned int
409 
413  unsigned int
415 
428  std::uint8_t
429  get_face_no(const unsigned int v = 0) const;
430 
440  unsigned int
442 
451  std::uint8_t
452  get_face_orientation(const unsigned int v = 0) const;
453 
463 
472  bool
474 
479  const std::array<unsigned int, n_lanes> &
480  get_cell_ids() const
481  {
482 // implemented inline to avoid compilation problems on Windows
483 #ifdef DEBUG
485 #endif
486  return cell_ids;
487  }
488 
493  const std::array<unsigned int, n_lanes> &
494  get_face_ids() const
495  {
496 // implemented inline to avoid compilation problems on Windows
497 #ifdef DEBUG
499 #endif
500  return face_ids;
501  }
502 
507  unsigned int
509  {
510 // implemented inline to avoid compilation problems on Windows
511 #ifdef DEBUG
513 #endif
514 
515  return cell;
516  }
517 
522  const std::array<unsigned int, n_lanes> &
524  {
525 // implemented inline to avoid compilation problems on Windows
526 #ifdef DEBUG
528 #endif
529 
530  if (!is_face || dof_access_index ==
532  return cell_ids;
533  else
534  return face_ids;
535  }
536 
544 
546 
551 
558  Number
560 
567  void
568  set_cell_data(AlignedVector<Number> &array, const Number &value) const;
569 
574  template <typename T>
575  std::array<T, Number::size()>
577  const AlignedVector<std::array<T, Number::size()>> &array) const;
578 
583  template <typename T>
584  void
585  set_cell_data(AlignedVector<std::array<T, Number::size()>> &array,
586  const std::array<T, Number::size()> & value) const;
587 
596  Number
598 
607  void
608  set_face_data(AlignedVector<Number> &array, const Number &value) const;
609 
614  template <typename T>
615  std::array<T, Number::size()>
617  const AlignedVector<std::array<T, Number::size()>> &array) const;
618 
623  template <typename T>
624  void
625  set_face_data(AlignedVector<std::array<T, Number::size()>> &array,
626  const std::array<T, Number::size()> & value) const;
627 
629 
635  {
637  const DoFInfo * dof_info;
639  unsigned int active_fe_index;
640  unsigned int active_quad_index;
642  };
643 
644 protected:
651  FEEvaluationData(const InitializationData &initialization_data,
652  const bool is_interior_face,
653  const unsigned int quad_no,
654  const unsigned int first_selected_component);
655 
661  const std::shared_ptr<
663  & mapping_data,
664  const unsigned int n_fe_components,
665  const unsigned int first_selected_component);
666 
674 
681 
689 
695  const unsigned int quad_no;
696 
701  const unsigned int n_fe_components;
702 
707  const unsigned int first_selected_component;
708 
712  const unsigned int active_fe_index;
713 
718  const unsigned int active_quad_index;
719 
726 
730  const unsigned int n_quadrature_points;
731 
737 
751 
756  const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>>
758 
765  const Number *J_value;
766 
771 
776 
781 
789 
800  Number *values_dofs;
801 
811  Number *values_quad;
812 
824  Number *gradients_quad;
825 
838 
848  Number *hessians_quad;
849 
855 
862 
869 
876 
883 
890 
897 
904 
909  unsigned int cell;
910 
918 
924 
933  std::array<std::uint8_t, n_lanes> face_numbers;
934 
943  std::array<std::uint8_t, n_lanes> face_orientations;
944 
952  unsigned int subface_index;
953 
961 
966  std::array<unsigned int, n_lanes> cell_ids;
967 
972  std::array<unsigned int, n_lanes> face_ids;
973 
979  std::shared_ptr<
982 
983  // Make FEEvaluation and FEEvaluationBase objects friends for access to
984  // protected member mapped_geometry.
985  template <int, int, typename, bool, typename>
986  friend class FEEvaluationBase;
987 
988  template <int, int, int, int, typename, typename>
989  friend class FEEvaluation;
990 };
991 
992 
997 template <int dim,
998  typename Number,
999  bool is_face,
1000  typename VectorizedArrayType = VectorizedArray<Number>>
1003 
1004 
1005 /*----------------------- Inline functions ----------------------------------*/
1006 
1007 #ifndef DOXYGEN
1008 
1009 template <int dim, typename Number, bool is_face>
1011  const ShapeInfoType &shape_info,
1012  const bool is_interior_face)
1013  : FEEvaluationData(
1014  InitializationData{&shape_info, nullptr, nullptr, 0, 0, nullptr},
1016  0,
1017  0)
1018 {}
1019 
1020 
1021 template <int dim, typename Number, bool is_face>
1023  const InitializationData &initialization_data,
1024  const bool is_interior_face,
1025  const unsigned int quad_no,
1026  const unsigned int first_selected_component)
1027  : data(initialization_data.shape_info)
1028  , dof_info(initialization_data.dof_info)
1029  , mapping_data(initialization_data.mapping_data)
1030  , quad_no(quad_no)
1031  , n_fe_components(dof_info != nullptr ? dof_info->start_components.back() : 0)
1033  , active_fe_index(initialization_data.active_fe_index)
1034  , active_quad_index(initialization_data.active_quad_index)
1035  , descriptor(initialization_data.descriptor)
1036  , n_quadrature_points(descriptor == nullptr ?
1037  (is_face ? data->n_q_points_face : data->n_q_points) :
1039  , quadrature_points(nullptr)
1040  , jacobian(nullptr)
1041  , jacobian_gradients(nullptr)
1042  , J_value(nullptr)
1043  , normal_vectors(nullptr)
1044  , normal_x_jacobian(nullptr)
1046  descriptor != nullptr ? descriptor->quadrature_weights.begin() : nullptr)
1047 # ifdef DEBUG
1048  , is_reinitialized(false)
1049  , dof_values_initialized(false)
1050  , values_quad_initialized(false)
1052  , hessians_quad_initialized(false)
1053  , values_quad_submitted(false)
1054  , gradients_quad_submitted(false)
1055 # endif
1058  , dof_access_index(
1059  is_face ?
1060  (is_interior_face ?
1061  internal::MatrixFreeFunctions::DoFInfo::dof_access_face_interior :
1062  internal::MatrixFreeFunctions::DoFInfo::dof_access_face_exterior) :
1063  internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
1064  , subface_index(0)
1065  , cell_type(internal::MatrixFreeFunctions::general)
1066 {}
1067 
1068 
1069 
1070 template <int dim, typename Number, bool is_face>
1072  const std::shared_ptr<
1074  & mapped_geometry,
1075  const unsigned int n_fe_components,
1076  const unsigned int first_selected_component)
1077  : data(nullptr)
1078  , dof_info(nullptr)
1079  , mapping_data(nullptr)
1085  , descriptor(nullptr)
1087  mapped_geometry->get_data_storage().descriptor[0].n_q_points)
1088  , quadrature_points(nullptr)
1089  , jacobian(nullptr)
1090  , jacobian_gradients(nullptr)
1091  , J_value(nullptr)
1092  , normal_vectors(nullptr)
1093  , normal_x_jacobian(nullptr)
1094  , quadrature_weights(nullptr)
1095  , cell(0)
1096  , cell_type(internal::MatrixFreeFunctions::general)
1097  , interior_face(true)
1098  , dof_access_index(internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
1100  , is_reinitialized(false)
1101 {
1102  mapping_data = &mapped_geometry->get_data_storage();
1103  jacobian = mapped_geometry->get_data_storage().jacobians[0].begin();
1104  J_value = mapped_geometry->get_data_storage().JxW_values.begin();
1106  mapped_geometry->get_data_storage().jacobian_gradients[0].begin();
1108  mapped_geometry->get_data_storage().quadrature_points.begin();
1109 }
1110 
1111 
1112 
1113 template <int dim, typename Number, bool is_face>
1116 {
1123 
1124  data = other.data;
1125  dof_info = other.dof_info;
1126  mapping_data = other.mapping_data;
1127  descriptor = other.descriptor;
1128  jacobian = nullptr;
1129  J_value = nullptr;
1130  normal_vectors = nullptr;
1131  normal_x_jacobian = nullptr;
1132  jacobian_gradients = nullptr;
1133  quadrature_points = nullptr;
1135 
1136 # ifdef DEBUG
1137  is_reinitialized = false;
1138  dof_values_initialized = false;
1139  values_quad_initialized = false;
1141  hessians_quad_initialized = false;
1142  values_quad_submitted = false;
1143  gradients_quad_submitted = false;
1144 # endif
1145 
1147  interior_face = other.is_interior_face();
1149  is_face ?
1150  (is_interior_face() ?
1154  face_numbers[0] = 0;
1155  face_orientations[0] = 0;
1156  subface_index = 0;
1158 
1159  return *this;
1160 }
1161 
1162 
1163 
1164 template <int dim, typename Number, bool is_face>
1165 inline void
1168  const unsigned int n_components)
1169 {
1171 
1172  const unsigned int tensor_dofs_per_component =
1173  Utilities::fixed_power<dim>(data->data.front().fe_degree + 1);
1174  const unsigned int dofs_per_component = data->dofs_per_component_on_cell;
1175 
1176  const unsigned int size_scratch_data =
1177  std::max(tensor_dofs_per_component + 1, dofs_per_component) * n_components *
1178  3 +
1179  2 * n_quadrature_points;
1180  const unsigned int size_data_arrays =
1182  (n_components * ((dim * (dim + 1)) / 2 + 2 * dim + 1) *
1184 
1185  const unsigned int allocated_size = size_scratch_data + size_data_arrays;
1186  scratch_data_array->resize_fast(allocated_size);
1187  scratch_data.reinit(scratch_data_array->begin() + size_data_arrays,
1188  size_scratch_data);
1189 
1190  // set the pointers to the correct position in the data array
1198  hessians_quad =
1200  n_components * (dofs_per_component + (2 * dim + 1) * n_quadrature_points);
1201 }
1202 
1203 
1204 
1205 template <int dim, typename Number, bool is_face>
1206 inline void
1209 {
1210  Assert(is_face == true,
1211  ExcMessage("Faces can only be set if the is_face template parameter "
1212  "is true"));
1213  face_numbers[0] =
1215  subface_index = is_interior_face() == true ?
1217  face.subface_index;
1218 
1219  // First check if interior or exterior cell has non-standard orientation
1220  // (i.e. the third bit is one or not). Then set zero if this cell has
1221  // standard-orientation else copy the first three bits
1222  // (which is equivalent to modulo 8). See also the documentation of
1223  // internal::MatrixFreeFunctions::FaceToCellTopology::face_orientation.
1224  face_orientations[0] = (is_interior_face() == (face.face_orientation >= 8)) ?
1225  (face.face_orientation % 8) :
1226  0;
1227 
1228  if (is_interior_face())
1229  cell_ids = face.cells_interior;
1230  else
1231  cell_ids = face.cells_exterior;
1232 }
1233 
1234 
1235 
1236 template <int dim, typename Number, bool is_face>
1239  const unsigned int q_point) const
1240 {
1242  Assert(normal_vectors != nullptr,
1244  "update_normal_vectors"));
1246  return normal_vectors[0];
1247  else
1248  return normal_vectors[q_point];
1249 }
1250 
1251 
1252 
1253 template <int dim, typename Number, bool is_face>
1254 inline DEAL_II_ALWAYS_INLINE Number
1255 FEEvaluationData<dim, Number, is_face>::JxW(const unsigned int q_point) const
1256 {
1258  Assert(J_value != nullptr,
1260  "update_values|update_gradients"));
1262  {
1264  return J_value[0] * quadrature_weights[q_point];
1265  }
1266  else
1267  return J_value[q_point];
1268 }
1269 
1270 
1271 
1272 template <int dim, typename Number, bool is_face>
1275  const unsigned int q) const
1276 {
1278  Assert(this->quadrature_points != nullptr,
1280  "update_quadrature_points"));
1281 
1282  // Cartesian/affine mesh: only first vertex of cell is stored, we must
1283  // compute it through the Jacobian (which is stored in non-inverted and
1284  // non-transposed form as index '1' in the jacobian field)
1285  if (is_face == false &&
1287  {
1288  Assert(this->jacobian != nullptr, ExcNotInitialized());
1290 
1291  const Tensor<2, dim, Number> &jac = this->jacobian[1];
1293  for (unsigned int d = 0; d < dim; ++d)
1294  point[d] += jac[d][d] * static_cast<typename Number::value_type>(
1295  this->descriptor->quadrature.point(q)[d]);
1296  else
1297  for (unsigned int d = 0; d < dim; ++d)
1298  for (unsigned int e = 0; e < dim; ++e)
1299  point[d] += jac[d][e] * static_cast<typename Number::value_type>(
1300  this->descriptor->quadrature.point(q)[e]);
1301  return point;
1302  }
1303  else
1304  return this->quadrature_points[q];
1305 }
1306 
1307 
1308 
1309 template <int dim, typename Number, bool is_face>
1312  const unsigned int q_point) const
1313 {
1315  Assert(jacobian != nullptr,
1317  "update_gradients"));
1319  return jacobian[0];
1320  else
1321  return jacobian[q_point];
1322 }
1323 
1324 
1325 
1326 template <int dim, typename Number, bool is_face>
1327 inline const Number *
1329 {
1330  return values_dofs;
1331 }
1332 
1333 
1334 
1335 template <int dim, typename Number, bool is_face>
1336 inline Number *
1338 {
1339 # ifdef DEBUG
1340  dof_values_initialized = true;
1341 # endif
1342  return values_dofs;
1343 }
1344 
1345 
1346 
1347 template <int dim, typename Number, bool is_face>
1348 inline const Number *
1350 {
1351 # ifdef DEBUG
1353 # endif
1354  return values_quad;
1355 }
1356 
1357 
1358 
1359 template <int dim, typename Number, bool is_face>
1360 inline Number *
1362 {
1363 # ifdef DEBUG
1364  values_quad_initialized = true;
1365  values_quad_submitted = true;
1366 # endif
1367  return values_quad;
1368 }
1369 
1370 
1371 
1372 template <int dim, typename Number, bool is_face>
1373 inline const Number *
1375 {
1376 # ifdef DEBUG
1378  ExcNotInitialized());
1379 # endif
1380  return gradients_quad;
1381 }
1382 
1383 
1384 
1385 template <int dim, typename Number, bool is_face>
1386 inline Number *
1388 {
1389 # ifdef DEBUG
1390  gradients_quad_submitted = true;
1392 # endif
1393  return gradients_quad;
1394 }
1395 
1396 
1397 
1398 template <int dim, typename Number, bool is_face>
1399 inline const Number *
1401 {
1402 # ifdef DEBUG
1404 # endif
1405  return hessians_quad;
1406 }
1407 
1408 
1409 
1410 template <int dim, typename Number, bool is_face>
1411 inline Number *
1413 {
1414 # ifdef DEBUG
1416 # endif
1417  return hessians_quad;
1418 }
1419 
1420 
1421 
1422 template <int dim, typename Number, bool is_face>
1423 inline unsigned int
1425 {
1426  Assert(mapping_data != nullptr, ExcInternalError());
1427 
1428  if (dof_info == nullptr)
1429  return 0;
1430  else
1431  {
1434  }
1435 }
1436 
1437 
1438 
1439 template <int dim, typename Number, bool is_face>
1442 {
1443 # ifdef DEBUG
1445 # endif
1446  return cell_type;
1447 }
1448 
1449 
1450 
1451 template <int dim, typename Number, bool is_face>
1454 {
1455  Assert(data != nullptr, ExcInternalError());
1456  return *data;
1457 }
1458 
1459 
1460 
1461 template <int dim, typename Number, bool is_face>
1464 {
1465  Assert(dof_info != nullptr,
1466  ExcMessage(
1467  "FEEvaluation was not initialized with a MatrixFree object!"));
1468  return *dof_info;
1469 }
1470 
1471 
1472 
1473 template <int dim, typename Number, bool is_face>
1474 inline const std::vector<unsigned int> &
1476 {
1477  return data->lexicographic_numbering;
1478 }
1479 
1480 
1481 
1482 template <int dim, typename Number, bool is_face>
1483 inline unsigned int
1485 {
1486  return quad_no;
1487 }
1488 
1489 
1490 
1491 template <int dim, typename Number, bool is_face>
1492 inline unsigned int
1494 {
1495  if (is_face && dof_access_index ==
1498  else
1499  return cell;
1500 }
1501 
1502 
1503 
1504 template <int dim, typename Number, bool is_face>
1505 inline unsigned int
1507 {
1508  return active_fe_index;
1509 }
1510 
1511 
1512 
1513 template <int dim, typename Number, bool is_face>
1514 inline unsigned int
1516 {
1517  return active_quad_index;
1518 }
1519 
1520 
1521 
1522 template <int dim, typename Number, bool is_face>
1523 inline unsigned int
1525 {
1526  return first_selected_component;
1527 }
1528 
1529 
1530 
1531 template <int dim, typename Number, bool is_face>
1532 inline ArrayView<Number>
1534 {
1535  return scratch_data;
1536 }
1537 
1538 
1539 
1540 template <int dim, typename Number, bool is_face>
1541 std::uint8_t
1542 FEEvaluationData<dim, Number, is_face>::get_face_no(const unsigned int v) const
1543 {
1544  Assert(is_face, ExcNotInitialized());
1545  Assert(v == 0 || (cell != numbers::invalid_unsigned_int &&
1546  dof_access_index ==
1548  is_interior_face() == false),
1549  ExcMessage("All face numbers can only be queried for ECL at exterior "
1550  "faces. Use get_face_no() in other cases."));
1551 
1552  return face_numbers[v];
1553 }
1554 
1555 
1556 
1557 template <int dim, typename Number, bool is_face>
1558 inline unsigned int
1560 {
1561  return subface_index;
1562 }
1563 
1564 
1565 
1566 template <int dim, typename Number, bool is_face>
1567 std::uint8_t
1569  const unsigned int v) const
1570 {
1571  Assert(is_face, ExcNotInitialized());
1572  Assert(v == 0 || (cell != numbers::invalid_unsigned_int &&
1573  dof_access_index ==
1575  is_interior_face() == false),
1576  ExcMessage("All face numbers can only be queried for ECL at exterior "
1577  "faces. Use get_face_no() in other cases."));
1578 
1579  return face_orientations[v];
1580 }
1581 
1582 
1583 
1584 template <int dim, typename Number, bool is_face>
1587 {
1588  return dof_access_index;
1589 }
1590 
1591 
1592 
1593 template <int dim, typename Number, bool is_face>
1594 inline bool
1596 {
1597  return interior_face;
1598 }
1599 
1600 
1601 
1602 template <int dim, typename Number, bool is_face>
1605 {
1606  return {0U, n_quadrature_points};
1607 }
1608 
1609 
1610 
1611 namespace internal
1612 {
1613  template <std::size_t N,
1614  typename VectorizedArrayType2,
1615  typename GlobalVectorType,
1616  typename FU>
1617  inline void
1618  process_cell_or_face_data(const std::array<unsigned int, N> indices,
1619  GlobalVectorType & array,
1620  VectorizedArrayType2 & out,
1621  const FU & fu)
1622  {
1623  for (unsigned int i = 0; i < N; ++i)
1624  if (indices[i] != numbers::invalid_unsigned_int)
1625  {
1626  AssertIndexRange(indices[i] / N, array.size());
1627  fu(out[i], array[indices[i] / N][indices[i] % N]);
1628  }
1629  }
1630 } // namespace internal
1631 
1632 
1633 
1634 template <int dim, typename Number, bool is_face>
1635 inline Number
1637  const AlignedVector<Number> &array) const
1638 {
1639  Number out = Number(1.);
1640  internal::process_cell_or_face_data(this->get_cell_ids(),
1641  array,
1642  out,
1643  [](auto &local, const auto &global) {
1644  local = global;
1645  });
1646  return out;
1647 }
1648 
1649 
1650 
1651 template <int dim, typename Number, bool is_face>
1652 inline void
1654  AlignedVector<Number> &array,
1655  const Number & in) const
1656 {
1657  internal::process_cell_or_face_data(this->get_cell_ids(),
1658  array,
1659  in,
1660  [](const auto &local, auto &global) {
1661  global = local;
1662  });
1663 }
1664 
1665 
1666 
1667 template <int dim, typename Number, bool is_face>
1668 template <typename T>
1669 inline std::array<T, Number::size()>
1671  const AlignedVector<std::array<T, Number::size()>> &array) const
1672 {
1673  std::array<T, Number::size()> out;
1674  internal::process_cell_or_face_data(this->get_cell_ids(),
1675  array,
1676  out,
1677  [](auto &local, const auto &global) {
1678  local = global;
1679  });
1680  return out;
1681 }
1682 
1683 
1684 
1685 template <int dim, typename Number, bool is_face>
1686 template <typename T>
1687 inline void
1689  AlignedVector<std::array<T, Number::size()>> &array,
1690  const std::array<T, Number::size()> & in) const
1691 {
1692  internal::process_cell_or_face_data(this->get_cell_ids(),
1693  array,
1694  in,
1695  [](const auto &local, auto &global) {
1696  global = local;
1697  });
1698 }
1699 
1700 
1701 
1702 template <int dim, typename Number, bool is_face>
1703 inline Number
1705  const AlignedVector<Number> &array) const
1706 {
1707  Number out = Number(1.);
1708  internal::process_cell_or_face_data(this->get_face_ids(),
1709  array,
1710  out,
1711  [](auto &local, const auto &global) {
1712  local = global;
1713  });
1714  return out;
1715 }
1716 
1717 
1718 
1719 template <int dim, typename Number, bool is_face>
1720 inline void
1722  AlignedVector<Number> &array,
1723  const Number & in) const
1724 {
1725  internal::process_cell_or_face_data(this->get_face_ids(),
1726  array,
1727  in,
1728  [](const auto &local, auto &global) {
1729  global = local;
1730  });
1731 }
1732 
1733 
1734 
1735 template <int dim, typename Number, bool is_face>
1736 template <typename T>
1737 inline std::array<T, Number::size()>
1739  const AlignedVector<std::array<T, Number::size()>> &array) const
1740 {
1741  std::array<T, Number::size()> out;
1742  internal::process_cell_or_face_data(this->get_face_ids(),
1743  array,
1744  out,
1745  [](auto &local, const auto &global) {
1746  local = global;
1747  });
1748  return out;
1749 }
1750 
1751 
1752 
1753 template <int dim, typename Number, bool is_face>
1754 template <typename T>
1755 inline void
1757  AlignedVector<std::array<T, Number::size()>> &array,
1758  const std::array<T, Number::size()> & in) const
1759 {
1760  internal::process_cell_or_face_data(this->get_face_ids(),
1761  array,
1762  in,
1763  [](const auto &local, auto &global) {
1764  global = local;
1765  });
1766 }
1767 
1768 
1769 
1770 #endif // ifndef DOXYGEN
1771 
1772 
1774 
1775 #endif
void resize_fast(const size_type new_size)
iterator begin()
size_type size() const
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
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
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