Reference documentation for deal.II version GIT relicensing-437-g81ec864850 2024-04-19 07:30: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\}}\)
Loading...
Searching...
No Matches
fe_evaluation_data.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2021 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15
16#ifndef dealii_matrix_free_fe_evaluation_data_h
17#define dealii_matrix_free_fe_evaluation_data_h
18
19
20#include <deal.II/base/config.h>
21
30#include <deal.II/base/tensor.h>
32
36
37
39
40
41
42namespace 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 "
52 << "you pass to MatrixFree::reinit() via MatrixFree::AdditionalData. "
53 << "Here, 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
59template <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>>
65class FEEvaluation;
66
67template <int dim,
68 int n_components_,
69 typename Number,
70 bool = false,
71 typename = VectorizedArray<Number>>
73
74namespace internal
75{
76 namespace MatrixFreeFunctions
77 {
78 template <int, typename>
79 class MappingDataOnTheFly;
80 }
81} // namespace internal
82
83
112template <int dim, typename Number, bool is_face>
114{
116 using MappingInfoStorageType = internal::MatrixFreeFunctions::
117 MappingInfoStorage<(is_face ? dim - 1 : dim), dim, Number>;
119
120public:
121 static constexpr unsigned int dimension = dim;
122
123 using NumberType = Number;
127
128 static constexpr unsigned int n_lanes =
130
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 normal_vector(const unsigned int q_point) const;
230
236 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT("Use normal_vector() instead.")
237 Tensor<1, dim, Number>
238 get_normal_vector(const unsigned int q_point) const;
239
254 const Number *
256
265 Number *
267
278 const Number *
280
291 Number *
293
305 const Number *
307
319 Number *
321
334 const Number *
336
349 Number *
351
367 unsigned int
369
376 internal::MatrixFreeFunctions::GeometryType
378
382 const ShapeInfoType &
384
388 const internal::MatrixFreeFunctions::DoFInfo &
390
396 const std::vector<unsigned int> &
398
402 unsigned int
404
408 unsigned int
410
415 unsigned int
417
422 unsigned int
424
428 unsigned int
430
443 std::uint8_t
444 get_face_no(const unsigned int v = 0) const;
445
455 unsigned int
457
466 std::uint8_t
467 get_face_orientation(const unsigned int v = 0) const;
468
476 internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex
478
487 bool
489
494 const std::array<unsigned int, n_lanes> &
496 {
497// implemented inline to avoid compilation problems on Windows
498#ifdef DEBUG
500#endif
501 return cell_ids;
502 }
503
508 const std::array<unsigned int, n_lanes> &
510 {
511// implemented inline to avoid compilation problems on Windows
512#ifdef DEBUG
514#endif
515 return face_ids;
516 }
517
522 unsigned int
524 {
525// implemented inline to avoid compilation problems on Windows
526#ifdef DEBUG
528#endif
529
530 return cell;
531 }
532
537 const std::array<unsigned int, n_lanes> &
539 {
540// implemented inline to avoid compilation problems on Windows
541#ifdef DEBUG
543#endif
544
545 if (!is_face || dof_access_index ==
547 return cell_ids;
548 else
549 return face_ids;
550 }
551
559
577 template <typename T>
578 T
579 read_cell_data(const AlignedVector<T> &array) const;
580
587 template <typename T>
588 void
589 set_cell_data(AlignedVector<T> &array, const T &value) const;
590
599 template <typename T>
600 T
601 read_face_data(const AlignedVector<T> &array) const;
602
611 template <typename T>
612 void
613 set_face_data(AlignedVector<T> &array, const T &value) const;
614
630
631protected:
638 FEEvaluationData(const InitializationData &initialization_data,
639 const bool is_interior_face,
640 const unsigned int quad_no,
641 const unsigned int first_selected_component);
642
648 const std::shared_ptr<
651 const unsigned int n_fe_components,
652 const unsigned int first_selected_component);
653
661
668
676
682 const unsigned int quad_no;
683
688 const unsigned int n_fe_components;
689
694 const unsigned int first_selected_component;
695
699 const unsigned int active_fe_index;
700
705 const unsigned int active_quad_index;
706
713
717 const unsigned int n_quadrature_points;
718
724
738
743 const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>>
745
750 const Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>>
752
759 const Number *J_value;
760
765
770
775
783
794 Number *values_dofs;
795
805 Number *values_quad;
806
820
833
846
857
863
870
877
884
891
898
905
912
917 unsigned int cell;
918
926
932
941 std::array<std::uint8_t, n_lanes> face_numbers;
942
951 std::array<std::uint8_t, n_lanes> face_orientations;
952
960 unsigned int subface_index;
961
969
974 std::array<unsigned int, n_lanes> cell_ids;
975
980 std::array<unsigned int, n_lanes> face_ids;
981
987 std::shared_ptr<
990
996
997 // Make FEEvaluation and FEEvaluationBase objects friends for access to
998 // protected member mapped_geometry.
999 template <int, int, typename, bool, typename>
1000 friend class FEEvaluationBase;
1001
1002 template <int, int, int, int, typename, typename>
1003 friend class FEEvaluation;
1004};
1005
1006
1011template <int dim,
1012 typename Number,
1013 bool is_face,
1014 typename VectorizedArrayType = VectorizedArray<Number>>
1017
1018
1019/*----------------------- Inline functions ----------------------------------*/
1020
1021#ifndef DOXYGEN
1022
1023template <int dim, typename Number, bool is_face>
1025 const ShapeInfoType &shape_info,
1026 const bool is_interior_face)
1028 InitializationData{&shape_info, nullptr, nullptr, 0, 0, nullptr},
1030 0,
1031 0)
1032{}
1033
1034
1035template <int dim, typename Number, bool is_face>
1037 const InitializationData &initialization_data,
1038 const bool is_interior_face,
1039 const unsigned int quad_no,
1040 const unsigned int first_selected_component)
1041 : data(initialization_data.shape_info)
1042 , dof_info(initialization_data.dof_info)
1043 , mapping_data(initialization_data.mapping_data)
1044 , quad_no(quad_no)
1045 , n_fe_components(dof_info != nullptr ? dof_info->start_components.back() : 0)
1047 , active_fe_index(initialization_data.active_fe_index)
1048 , active_quad_index(initialization_data.active_quad_index)
1049 , descriptor(initialization_data.descriptor)
1050 , n_quadrature_points(descriptor == nullptr ?
1051 (is_face ? data->n_q_points_face : data->n_q_points) :
1053 , quadrature_points(nullptr)
1054 , jacobian(nullptr)
1055 , jacobian_gradients(nullptr)
1057 , J_value(nullptr)
1058 , normal_vectors(nullptr)
1059 , normal_x_jacobian(nullptr)
1061 descriptor != nullptr ? descriptor->quadrature_weights.begin() : nullptr)
1062# ifdef DEBUG
1063 , is_reinitialized(false)
1064 , dof_values_initialized(false)
1068 , values_quad_submitted(false)
1070# endif
1071 , cell(numbers::invalid_unsigned_int)
1074 is_face ?
1076 internal::MatrixFreeFunctions::DoFInfo::dof_access_face_interior :
1077 internal::MatrixFreeFunctions::DoFInfo::dof_access_face_exterior) :
1078 internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
1079 , subface_index(0)
1080 , cell_type(internal::MatrixFreeFunctions::general)
1082{}
1083
1084
1085
1086template <int dim, typename Number, bool is_face>
1088 const std::shared_ptr<
1091 const unsigned int n_fe_components,
1092 const unsigned int first_selected_component)
1093 : data(nullptr)
1094 , dof_info(nullptr)
1095 , mapping_data(nullptr)
1096 , quad_no(numbers::invalid_unsigned_int)
1099 , active_fe_index(numbers::invalid_unsigned_int)
1100 , active_quad_index(numbers::invalid_unsigned_int)
1101 , descriptor(nullptr)
1103 mapped_geometry->get_data_storage().descriptor[0].n_q_points)
1104 , quadrature_points(nullptr)
1105 , jacobian(nullptr)
1106 , jacobian_gradients(nullptr)
1108 , J_value(nullptr)
1109 , normal_vectors(nullptr)
1110 , normal_x_jacobian(nullptr)
1111 , quadrature_weights(nullptr)
1112 , cell(0)
1113 , cell_type(internal::MatrixFreeFunctions::general)
1114 , interior_face(true)
1115 , dof_access_index(internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
1117 , is_reinitialized(false)
1119{
1120 mapping_data = &mapped_geometry->get_data_storage();
1121 jacobian = mapped_geometry->get_data_storage().jacobians[0].begin();
1122 J_value = mapped_geometry->get_data_storage().JxW_values.begin();
1124 mapped_geometry->get_data_storage().jacobian_gradients[0].begin();
1126 .jacobian_gradients_non_inverse[0]
1127 .begin();
1129 mapped_geometry->get_data_storage().quadrature_points.begin();
1130}
1131
1132
1133
1134template <int dim, typename Number, bool is_face>
1137{
1144
1145 data = other.data;
1146 dof_info = other.dof_info;
1147 mapping_data = other.mapping_data;
1148 descriptor = other.descriptor;
1149 jacobian = nullptr;
1150 J_value = nullptr;
1151 normal_vectors = nullptr;
1152 normal_x_jacobian = nullptr;
1153 jacobian_gradients = nullptr;
1155 quadrature_points = nullptr;
1157
1158# ifdef DEBUG
1159 is_reinitialized = false;
1160 dof_values_initialized = false;
1164 values_quad_submitted = false;
1166# endif
1167
1171 is_face ?
1172 (is_interior_face() ?
1175 internal::MatrixFreeFunctions::DoFInfo::dof_access_cell;
1176 face_numbers[0] = 0;
1177 face_orientations[0] = 0;
1178 subface_index = 0;
1181
1182 return *this;
1183}
1184
1185
1186
1187template <int dim, typename Number, bool is_face>
1188inline void
1191 const unsigned int n_components)
1192{
1194
1195 const unsigned int tensor_dofs_per_component =
1196 Utilities::fixed_power<dim>(data->data.front().fe_degree + 1);
1198
1199 const unsigned int size_scratch_data =
1200 std::max(tensor_dofs_per_component + 1, dofs_per_component) * n_components *
1201 3 +
1203 const unsigned int size_data_arrays =
1205 (n_components * ((dim * (dim + 1)) / 2 + 2 * dim + 2) *
1207
1208 // include 12 extra fields to insert some padding between values, gradients
1209 // and hessians, which helps to reduce the probability of cache conflicts
1210 const unsigned int allocated_size = size_scratch_data + size_data_arrays + 12;
1211# ifdef DEBUG
1213 scratch_data_array->resize(allocated_size,
1214 Number(numbers::signaling_nan<ScalarNumber>()));
1215# else
1216 scratch_data_array->resize_fast(allocated_size);
1217# endif
1218 scratch_data.reinit(scratch_data_array->begin() + size_data_arrays + 12,
1219 size_scratch_data);
1220
1221 // set the pointers to the correct position in the data array
1223 values_quad =
1226 scratch_data_array->begin() + 6 +
1229 scratch_data_array->begin() + 8 +
1232 scratch_data_array->begin() + 12 +
1235 scratch_data_array->begin() + 12 +
1237}
1238
1239
1240
1241template <int dim, typename Number, bool is_face>
1242inline void
1245{
1246 Assert(is_face == true,
1247 ExcMessage("Faces can only be set if the is_face template parameter "
1248 "is true"));
1249 face_numbers[0] =
1251 subface_index = is_interior_face() == true ?
1253 face.subface_index;
1254
1255 // First check if interior or exterior cell has non-standard orientation
1256 // (i.e. the third bit is one or not). Then set zero if this cell has
1257 // standard-orientation else copy the first three bits
1258 // (which is equivalent to modulo 8). See also the documentation of
1259 // internal::MatrixFreeFunctions::FaceToCellTopology::face_orientation.
1260 face_orientations[0] = (is_interior_face() == (face.face_orientation >= 8)) ?
1261 (face.face_orientation % 8) :
1262 0;
1263
1264 if (is_interior_face())
1265 cell_ids = face.cells_interior;
1266 else
1267 cell_ids = face.cells_exterior;
1268}
1269
1270
1271
1272template <int dim, typename Number, bool is_face>
1275 const unsigned int q_point) const
1276{
1278 Assert(normal_vectors != nullptr,
1280 "update_normal_vectors"));
1282 return normal_vectors[0];
1283 else
1284 return normal_vectors[q_point];
1285}
1286
1287
1288
1289// This function is deprecated.
1290template <int dim, typename Number, bool is_face>
1293 const unsigned int q_point) const
1294{
1295 return normal_vector(q_point);
1296}
1297
1298
1299
1300template <int dim, typename Number, bool is_face>
1301inline DEAL_II_ALWAYS_INLINE Number
1302FEEvaluationData<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
1319template <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
1356template <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
1373template <int dim, typename Number, bool is_face>
1374inline const Number *
1376{
1377 return values_dofs;
1378}
1379
1380
1381
1382template <int dim, typename Number, bool is_face>
1383inline Number *
1385{
1386# ifdef DEBUG
1388# endif
1389 return values_dofs;
1390}
1391
1392
1393
1394template <int dim, typename Number, bool is_face>
1395inline const Number *
1397{
1398# ifdef DEBUG
1400# endif
1401 return values_quad;
1402}
1403
1404
1405
1406template <int dim, typename Number, bool is_face>
1407inline Number *
1409{
1410# ifdef DEBUG
1412 values_quad_submitted = true;
1413# endif
1414 return values_quad;
1415}
1416
1417
1418
1419template <int dim, typename Number, bool is_face>
1420inline const Number *
1422{
1423# ifdef DEBUG
1426# endif
1427 return gradients_quad;
1428}
1429
1430
1431
1432template <int dim, typename Number, bool is_face>
1433inline Number *
1435{
1436# ifdef DEBUG
1439# endif
1440 return gradients_quad;
1441}
1442
1443
1444
1445template <int dim, typename Number, bool is_face>
1446inline const Number *
1448{
1449# ifdef DEBUG
1451# endif
1452 return hessians_quad;
1453}
1454
1455
1456
1457template <int dim, typename Number, bool is_face>
1458inline Number *
1460{
1461# ifdef DEBUG
1463# endif
1464 return hessians_quad;
1465}
1466
1467
1468
1469template <int dim, typename Number, bool is_face>
1470inline 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
1486template <int dim, typename Number, bool is_face>
1489{
1490# ifdef DEBUG
1492# endif
1493 return cell_type;
1494}
1495
1496
1497
1498template <int dim, typename Number, bool is_face>
1501{
1502 Assert(data != nullptr, ExcInternalError());
1503 return *data;
1504}
1505
1506
1507
1508template <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
1520template <int dim, typename Number, bool is_face>
1521inline const std::vector<unsigned int> &
1523{
1525}
1526
1527
1528
1529template <int dim, typename Number, bool is_face>
1530inline unsigned int
1532{
1533 return quad_no;
1534}
1535
1536
1537
1538template <int dim, typename Number, bool is_face>
1539inline unsigned int
1541{
1542 if (is_face && dof_access_index ==
1545 else
1546 return cell;
1547}
1548
1549
1550
1551template <int dim, typename Number, bool is_face>
1552inline unsigned int
1554{
1555 return active_fe_index;
1556}
1557
1558
1559
1560template <int dim, typename Number, bool is_face>
1561inline unsigned int
1563{
1564 return active_quad_index;
1565}
1566
1567
1568
1569template <int dim, typename Number, bool is_face>
1570inline unsigned int
1572{
1574}
1575
1576
1577
1578template <int dim, typename Number, bool is_face>
1579inline ArrayView<Number>
1581{
1582 return scratch_data;
1583}
1584
1585
1586
1587template <int dim, typename Number, bool is_face>
1588inline std::uint8_t
1590{
1591 Assert(is_face, ExcNotInitialized());
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
1604template <int dim, typename Number, bool is_face>
1605inline unsigned int
1607{
1608 return subface_index;
1609}
1610
1611
1612
1613template <int dim, typename Number, bool is_face>
1614inline std::uint8_t
1616 const unsigned int v) const
1617{
1618 Assert(is_face, ExcNotInitialized());
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
1631template <int dim, typename Number, bool is_face>
1634{
1635 return dof_access_index;
1636}
1637
1638
1639
1640template <int dim, typename Number, bool is_face>
1641inline bool
1643{
1644 return interior_face;
1645}
1646
1647
1648
1649template <int dim, typename Number, bool is_face>
1652{
1653 return {0U, n_quadrature_points};
1654}
1655
1656
1657
1658namespace internal
1659{
1660 template <std::size_t N,
1661 typename VectorOfArrayType,
1662 typename ArrayType,
1663 typename FU>
1664 inline void
1665 process_cell_or_face_data(const std::array<unsigned int, N> indices,
1666 VectorOfArrayType &array,
1667 ArrayType &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
1678 template <std::size_t N, typename VectorOfArrayType, typename ArrayType>
1679 inline void
1680 set_valid_element_to_array(const std::array<unsigned int, N> indices,
1681 const VectorOfArrayType &array,
1682 ArrayType &out)
1683 {
1684 AssertDimension(indices.size(), out.size());
1685 AssertDimension(indices.size(), array[0].size());
1686 // set all entries in array to a valid element
1687 int index = 0;
1688 for (; index < N; ++index)
1689 if (indices[index] != numbers::invalid_unsigned_int)
1690 break;
1691 for (unsigned int i = 0; i < N; ++i)
1692 out[i] = array[indices[index] / N][indices[index] % N];
1693 }
1694} // namespace internal
1695
1696
1697
1698template <int dim, typename Number, bool is_face>
1699template <typename T>
1700inline T
1702 const AlignedVector<T> &array) const
1703{
1704 T out;
1705 internal::set_valid_element_to_array(this->get_cell_ids(), array, out);
1706 internal::process_cell_or_face_data(this->get_cell_ids(),
1707 array,
1708 out,
1709 [](auto &local, const auto &global) {
1710 local = global;
1711 });
1712 return out;
1713}
1714
1715
1716
1717template <int dim, typename Number, bool is_face>
1718template <typename T>
1719inline void
1721 const T &in) const
1722{
1723 internal::process_cell_or_face_data(this->get_cell_ids(),
1724 array,
1725 in,
1726 [](const auto &local, auto &global) {
1727 global = local;
1728 });
1729}
1730
1731
1732
1733template <int dim, typename Number, bool is_face>
1734template <typename T>
1735inline T
1737 const AlignedVector<T> &array) const
1738{
1739 T out;
1740 internal::set_valid_element_to_array(this->get_cell_ids(), array, out);
1741 internal::process_cell_or_face_data(this->get_face_ids(),
1742 array,
1743 out,
1744 [](auto &local, const auto &global) {
1745 local = global;
1746 });
1747 return out;
1748}
1749
1750
1751
1752template <int dim, typename Number, bool is_face>
1753template <typename T>
1754inline void
1756 const T &in) const
1757{
1758 internal::process_cell_or_face_data(this->get_face_ids(),
1759 array,
1760 in,
1761 [](const auto &local, auto &global) {
1762 global = local;
1763 });
1764}
1765
1766
1767
1768#endif // ifndef DOXYGEN
1769
1770
1772
1773#endif
void resize_fast(const size_type new_size)
iterator begin()
size_type size() const
void resize(const size_type new_size)
ArrayView(const std::initializer_list< std::remove_cv_t< value_type > > &initializer_list) DEAL_II_CXX20_REQUIRES(std void reinit(value_type *starting_element, const std::size_t n_elements)
Definition array_view.h:344
const unsigned int quad_no
const Number * J_value
internal::MatrixFreeFunctions::GeometryType get_cell_type() const
const Tensor< 2, dim, Number > * jacobian
const MappingInfoStorageType::QuadratureDescriptor * descriptor
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
const ShapeInfoType & get_shape_info() const
void set_face_data(AlignedVector< T > &array, const T &value) const
void reinit_face(const internal::MatrixFreeFunctions::FaceToCellTopology< n_lanes > &face)
Number JxW(const unsigned int q_point) const
const std::array< unsigned int, n_lanes > & get_face_ids() const
unsigned int get_first_selected_component() const
const unsigned int n_fe_components
const ShapeInfoType * data
T read_face_data(const AlignedVector< T > &array) const
Number * gradients_from_hessians_quad
std::shared_ptr< internal::MatrixFreeFunctions::MappingDataOnTheFly< dim, Number > > mapped_geometry
unsigned int get_active_quadrature_index() const
const internal::MatrixFreeFunctions::DoFInfo & get_dof_info() const
const std::array< unsigned int, n_lanes > & get_cell_ids() const
unsigned int get_mapping_data_index_offset() const
void set_data_pointers(AlignedVector< Number > *scratch_data, const unsigned int n_components)
Tensor< 2, dim, Number > inverse_jacobian(const unsigned int q_point) const
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
internal::MatrixFreeFunctions::GeometryType cell_type
const Number * begin_gradients() const
const std::array< unsigned int, n_lanes > & get_cell_or_face_ids() const
Tensor< 1, dim, Number > normal_vector(const unsigned int q_point) const
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
Tensor< 1, dim, Number > get_normal_vector(const unsigned int q_point) const
const unsigned int active_quad_index
unsigned int get_active_fe_index() const
Number * values_from_gradients_quad
void set_cell_data(AlignedVector< T > &array, const T &value) 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
const Tensor< 1, dim, Number > * normal_x_jacobian
bool is_interior_face() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > quadrature_point_indices() const
const std::vector< unsigned int > & get_internal_dof_numbering() const
FEEvaluationData & operator=(const FEEvaluationData &other)
virtual ~FEEvaluationData()=default
ArrayView< Number > get_scratch_data() const
T read_cell_data(const AlignedVector< T > &array) 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)
FEEvaluationData(const ShapeInfoType &shape_info, const bool is_interior_face=true)
Point< dim, Number > quadrature_point(const unsigned int q) const
const Number * begin_values() const
typename internal::VectorizedArrayTrait< Number >::value_type ScalarNumber
std::array< unsigned int, n_lanes > face_ids
FEEvaluationData(const InitializationData &initialization_data, const bool is_interior_face, const unsigned int quad_no, const unsigned int first_selected_component)
const Number * begin_dof_values() const
std::uint8_t get_face_orientation(const unsigned int v=0) const
const Number * begin_hessians() const
unsigned int get_cell_or_face_batch_id() const
const unsigned int first_selected_component
const unsigned int dofs_per_component
const unsigned int n_q_points
static constexpr unsigned int n_components
Definition point.h:111
const Point< dim > & point(const unsigned int i) const
#define DEAL_II_ALWAYS_INLINE
Definition config.h:109
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
#define DeclException0(Exception0)
Definition exceptions.h:471
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcAccessToUninitializedField()
static ::ExceptionBase & ExcNotInitialized()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:516
static ::ExceptionBase & ExcMatrixFreeAccessToUninitializedMappingField(std::string arg1)
static ::ExceptionBase & ExcMessage(std::string arg1)
@ general
No special properties.
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition utilities.cc:191
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static const unsigned int invalid_unsigned_int
Definition types.h:220
const InputIterator OutputIterator out
Definition parallel.h:167
const Iterator & begin
Definition parallel.h:609
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:45
STL namespace.
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
const MappingInfoStorageType * mapping_data
const MappingInfoStorageType::QuadratureDescriptor * descriptor
std::array< unsigned int, vectorization_width > cells_interior
Definition face_info.h:60
std::array< unsigned int, vectorization_width > cells_exterior
Definition face_info.h:75
std::vector< UnivariateShapeData< Number > > data
Definition shape_info.h:485
std::vector< unsigned int > lexicographic_numbering
Definition shape_info.h:479
static constexpr std::size_t width()