deal.II version GIT relicensing-2394-g350a75cc75 2025-01-21 18:30:00+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
reference_cell.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2020 - 2024 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#ifndef dealii_tria_reference_cell_h
16#define dealii_tria_reference_cell_h
17
18#include <deal.II/base/config.h>
19
20#include <deal.II/base/array_view.h>
21#include <deal.II/base/geometry_info.h>
22#include <deal.II/base/ndarray.h>
23#include <deal.II/base/point.h>
24#include <deal.II/base/tensor.h>
25#include <deal.II/base/utilities.h>
26
27#include <deal.II/grid/tria_orientation.h>
28
29#include <boost/container/small_vector.hpp>
30
31#include <iosfwd>
32#include <string>
33
34DEAL_II_NAMESPACE_OPEN
35
36// Forward declarations
37#ifndef DOXYGEN
38template <int dim, int spacedim>
39class Mapping;
40
41template <int dim>
42class Quadrature;
43
44class ReferenceCell;
45#endif
46
47
48namespace internal
49{
62 constexpr ReferenceCell
63 make_reference_cell_from_int(const std::uint8_t kind);
64} // namespace internal
65
66
67
74enum class IsotropicRefinementChoice : std::uint8_t
75{
79 isotropic_refinement = 0,
83 cut_tet_68 = 1,
87 cut_tet_57 = 2,
91 cut_tet_49 = 3,
92};
93
94
95
126class ReferenceCell
127{
128public:
136 static ReferenceCell
137 n_vertices_to_type(const int dim, const unsigned int n_vertices);
138
147 constexpr ReferenceCell();
148
159 bool
160 is_hyper_cube() const;
161
165 bool
166 is_simplex() const;
167
172 unsigned int
173 get_dimension() const;
174
192 template <int dim>
193 double
194 d_linear_shape_function(const Point<dim> &xi, const unsigned int i) const;
195
204 template <int dim>
205 Tensor<1, dim>
206 d_linear_shape_function_gradient(const Point<dim> &xi,
207 const unsigned int i) const;
208
218 template <int dim, int spacedim = dim>
219 std::unique_ptr<Mapping<dim, spacedim>>
220 get_default_mapping(const unsigned int degree) const;
221
233 template <int dim, int spacedim = dim>
234 const Mapping<dim, spacedim> &
235 get_default_linear_mapping() const;
236
249 template <int dim>
250 Quadrature<dim>
251 get_gauss_type_quadrature(const unsigned n_points_1d) const;
252
266 template <int dim>
267 Quadrature<dim>
268 get_midpoint_quadrature() const;
269
287 template <int dim>
288 const Quadrature<dim> &
289 get_nodal_type_quadrature() const;
290
305 unsigned int
306 n_vertices() const;
307
312 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
313 vertex_indices() const;
314
326 template <int dim>
327 Point<dim>
328 vertex(const unsigned int v) const;
329
335 unsigned int
336 n_lines() const;
337
342 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
343 line_indices() const;
344
350 unsigned int
351 n_faces() const;
352
357 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
358 face_indices() const;
359
363 unsigned int
364 n_isotropic_refinement_choices() const;
365
369 IsotropicRefinementChoice
370 get_isotropic_refinement_choice(const unsigned int ref_choice) const;
371
383 unsigned int
384 n_isotropic_children() const;
385
394 template <int dim>
395 unsigned int
396 n_children(const RefinementCase<dim> ref_case =
397 RefinementCase<dim>::isotropic_refinement) const;
398
403 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
404 isotropic_child_indices() const;
405
418 ReferenceCell
419 face_reference_cell(const unsigned int face_no) const;
420
434 static constexpr types::geometric_orientation
435 default_combined_face_orientation();
436
447 static constexpr types::geometric_orientation
448 reversed_combined_line_orientation();
449
491 unsigned int
492 child_cell_on_face(const unsigned int face,
493 const unsigned int subface,
494 const types::geometric_orientation face_orientation =
495 default_combined_face_orientation()) const;
496
506 std::array<unsigned int, 2>
507 standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const;
508
518 std::array<unsigned int, 2>
519 standard_line_to_face_and_line_index(const unsigned int line) const;
520
529 unsigned int
530 line_to_cell_vertices(const unsigned int line,
531 const unsigned int vertex) const;
532
536 unsigned int
537 face_to_cell_lines(const unsigned int face,
538 const unsigned int line,
539 const types::geometric_orientation face_orientation) const;
540
544 unsigned int
545 face_to_cell_vertices(
546 const unsigned int face,
547 const unsigned int vertex,
548 const types::geometric_orientation face_orientation) const;
549
573 template <int dim>
574 Point<dim>
575 face_vertex_location(const unsigned int face,
576 const unsigned int vertex) const;
577
581 unsigned int
582 standard_to_real_face_vertex(
583 const unsigned int vertex,
584 const unsigned int face,
585 const types::geometric_orientation face_orientation) const;
586
590 unsigned int
591 standard_to_real_face_line(
592 const unsigned int line,
593 const unsigned int face,
594 const types::geometric_orientation face_orientation) const;
595
606 bool
607 standard_vs_true_line_orientation(
608 const unsigned int line,
609 const unsigned int face,
610 const types::geometric_orientation face_orientation,
611 const bool line_orientation) const;
612
636 double
637 volume() const;
638
655 template <int dim>
656 Point<dim>
657 barycenter() const;
658
682 template <int dim>
683 bool
684 contains_point(const Point<dim> &p, const double tolerance = 0) const;
685
694 template <int dim>
695 Point<dim>
696 closest_point(const Point<dim> &p) const;
697
709 template <int dim>
710 Tensor<1, dim>
711 unit_tangential_vectors(const unsigned int face_no,
712 const unsigned int i) const;
713
721 template <int dim>
722 Tensor<1, dim>
723 unit_normal_vectors(const unsigned int face_no) const;
724
730 unsigned int
731 n_face_orientations(const unsigned int face_no) const;
732
749 template <typename T, std::size_t N>
750 DEAL_II_DEPRECATED types::geometric_orientation
751 compute_orientation(const std::array<T, N> &vertices_0,
752 const std::array<T, N> &vertices_1) const;
753
779 template <typename T>
780 types::geometric_orientation
781 get_combined_orientation(const ArrayView<const T> &vertices_0,
782 const ArrayView<const T> &vertices_1) const;
783
784
798 template <typename T, std::size_t N>
799 DEAL_II_DEPRECATED std::array<T, N>
800 permute_according_orientation(const std::array<T, N> &vertices,
801 const unsigned int orientation) const;
802
814 template <typename T>
815 boost::container::small_vector<T, 8>
816 permute_by_combined_orientation(
817 const ArrayView<const T> &vertices,
818 const types::geometric_orientation orientation) const;
819
825 types::geometric_orientation
826 get_inverse_combined_orientation(
827 const types::geometric_orientation orientation) const;
828
832 ArrayView<const unsigned int>
833 faces_for_given_vertex(const unsigned int vertex_index) const;
834
839 constexpr ::ndarray<unsigned int, 12, 4>
840 new_isotropic_child_face_lines(const unsigned int refinement_choice) const;
841
846 constexpr ::ndarray<unsigned int, 12, 4, 2>
847 new_isotropic_child_face_line_vertices(
848 const unsigned int refinement_choice) const;
849
854 constexpr ::ndarray<unsigned int, 8, 6>
855 new_isotropic_child_cell_faces(const unsigned int refinement_choice) const;
856
861 constexpr ::ndarray<unsigned int, 8, 4>
862 new_isotropic_child_cell_vertices(const unsigned int refinement_choice) const;
863
876 unsigned int
877 exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const;
878
882 unsigned int
883 exodusii_face_to_deal_face(const unsigned int face_n) const;
884
888 unsigned int
889 unv_vertex_to_deal_vertex(const unsigned int vertex_n) const;
890
894 unsigned int
895 vtk_linear_type() const;
896
901 unsigned int
902 vtk_quadratic_type() const;
903
908 unsigned int
909 vtk_lagrange_type() const;
910
927 template <int dim>
928 unsigned int
929 vtk_lexicographic_to_node_index(
930 const std::array<unsigned, dim> &node_indices,
931 const std::array<unsigned, dim> &nodes_per_direction,
932 const bool legacy_format) const;
933
937 unsigned int
938 vtk_vertex_to_deal_vertex(const unsigned int vertex_index) const;
939
943 unsigned int
944 gmsh_element_type() const;
945
959 std::string
960 to_string() const;
961
965 constexpr operator std::uint8_t() const;
966
970 constexpr bool
971 operator==(const ReferenceCell &type) const;
972
976 constexpr bool
977 operator!=(const ReferenceCell &type) const;
978
984 template <class Archive>
985 void
986 serialize(Archive &archive, const unsigned int /*version*/);
987
991 static constexpr std::size_t
992 memory_consumption();
993
998private:
1002 std::uint8_t kind;
1003
1010 constexpr ReferenceCell(const std::uint8_t kind);
1011
1017 static constexpr ndarray<unsigned int, 1, 1> vertex_vertex_permutations = {
1018 {{{0}}}};
1019
1023 static constexpr ndarray<unsigned int, 2, 2> line_vertex_permutations = {
1024 {{{1, 0}}, {{0, 1}}}};
1025
1029 static constexpr ndarray<unsigned int, 6, 3> triangle_vertex_permutations = {
1030 {{{0, 2, 1}},
1031 {{0, 1, 2}},
1032 {{2, 1, 0}},
1033 {{2, 0, 1}},
1034 {{1, 0, 2}},
1035 {{1, 2, 0}}}};
1036
1040 static constexpr ndarray<unsigned int, 8, 4>
1041 quadrilateral_vertex_permutations = {{{{0, 2, 1, 3}},
1042 {{0, 1, 2, 3}},
1043 {{2, 3, 0, 1}},
1044 {{2, 0, 3, 1}},
1045 {{3, 1, 2, 0}},
1046 {{3, 2, 1, 0}},
1047 {{1, 0, 3, 2}},
1048 {{1, 3, 0, 2}}}};
1049
1054 friend constexpr ReferenceCell
1055 internal::make_reference_cell_from_int(const std::uint8_t);
1056
1057 friend std::ostream &
1058 operator<<(std::ostream &out, const ReferenceCell &reference_cell);
1059
1060 friend std::istream &
1061 operator>>(std::istream &in, ReferenceCell &reference_cell);
1062};
1063
1064
1072std::ostream &
1073operator<<(std::ostream &out, const ReferenceCell &reference_cell);
1074
1082std::istream &
1083operator>>(std::istream &in, ReferenceCell &reference_cell);
1084
1085
1086inline constexpr ReferenceCell::ReferenceCell(const std::uint8_t kind)
1087 : kind(kind)
1088{}
1089
1090
1091
1092inline constexpr ReferenceCell::operator std::uint8_t() const
1093{
1094 return kind;
1095}
1096
1097
1098
1099inline constexpr bool
1100ReferenceCell::operator==(const ReferenceCell &type) const
1101{
1102 return kind == type.kind;
1103}
1104
1105
1106
1107inline constexpr bool
1108ReferenceCell::operator!=(const ReferenceCell &type) const
1109{
1110 return kind != type.kind;
1111}
1112
1113
1114
1115namespace internal
1116{
1117 inline constexpr ReferenceCell
1118 make_reference_cell_from_int(const std::uint8_t kind)
1119 {
1120#ifndef DEAL_II_CXX14_CONSTEXPR_BUG
1121 // Make sure these are the only indices from which objects can be
1122 // created.
1123 Assert((kind == std::numeric_limits<std::uint8_t>::max()) || (kind < 8),
1124 ExcInternalError());
1125#endif
1126
1127 // Call the private constructor, which we can from here because this
1128 // function is a 'friend'.
1129 return {kind};
1130 }
1131} // namespace internal
1132
1133
1134
1142namespace ReferenceCells
1143{
1144 constexpr const ReferenceCell Vertex =
1145 internal::make_reference_cell_from_int(0);
1146 constexpr const ReferenceCell Line =
1147 internal::make_reference_cell_from_int(1);
1148 constexpr const ReferenceCell Triangle =
1149 internal::make_reference_cell_from_int(2);
1150 constexpr const ReferenceCell Quadrilateral =
1151 internal::make_reference_cell_from_int(3);
1152 constexpr const ReferenceCell Tetrahedron =
1153 internal::make_reference_cell_from_int(4);
1154 constexpr const ReferenceCell Pyramid =
1155 internal::make_reference_cell_from_int(5);
1156 constexpr const ReferenceCell Wedge =
1157 internal::make_reference_cell_from_int(6);
1158 constexpr const ReferenceCell Hexahedron =
1159 internal::make_reference_cell_from_int(7);
1160 constexpr const ReferenceCell Invalid =
1161 internal::make_reference_cell_from_int(
1162 std::numeric_limits<std::uint8_t>::max());
1163
1169 template <int dim>
1170 constexpr const ReferenceCell &
1171 get_simplex();
1172
1178 template <int dim>
1179 constexpr const ReferenceCell &
1180 get_hypercube();
1181
1189 template <int structdim>
1190 constexpr unsigned int
1191 max_n_vertices();
1192
1200 template <int structdim>
1201 constexpr unsigned int
1202 max_n_lines();
1203
1219 template <int structdim>
1220 constexpr unsigned int
1221 max_n_faces();
1222} // namespace ReferenceCells
1223
1224
1225
1226inline constexpr ReferenceCell::ReferenceCell()
1227 : ReferenceCell(ReferenceCells::Invalid)
1228{}
1229
1230
1231
1232template <class Archive>
1233inline void
1234ReferenceCell::serialize(Archive &archive, const unsigned int /*version*/)
1235{
1236 archive &kind;
1237}
1238
1239
1240
1241inline constexpr std::size_t
1246
1247
1248
1249inline ArrayView<const unsigned int>
1250ReferenceCell::faces_for_given_vertex(const unsigned int vertex) const
1251{
1252 AssertIndexRange(vertex, n_vertices());
1253 switch (this->kind)
1254 {
1255 case ReferenceCells::Line:
1256 return {&GeometryInfo<2>::vertex_to_face[vertex][0], 1};
1257 case ReferenceCells::Quadrilateral:
1258 return {&GeometryInfo<2>::vertex_to_face[vertex][0], 2};
1259 case ReferenceCells::Triangle:
1260 {
1261 static constexpr ndarray<unsigned int, 3, 2> table = {
1262 {{{0, 2}}, {{0, 1}}, {{1, 2}}}};
1263 return table[vertex];
1264 }
1265 case ReferenceCells::Tetrahedron:
1266 {
1267 static constexpr ndarray<unsigned int, 4, 3> table = {
1268 {{{0, 1, 2}}, {{0, 1, 3}}, {{0, 2, 3}}, {{1, 2, 3}}}};
1269
1270 return table[vertex];
1271 }
1272 case ReferenceCells::Pyramid:
1273 {
1274 static constexpr unsigned int X = numbers::invalid_unsigned_int;
1275 static constexpr ndarray<unsigned int, 5, 4> table = {
1276 {{{0, 1, 3, X}},
1277 {{0, 2, 3, X}},
1278 {{0, 1, 4, X}},
1279 {{0, 2, 4, X}},
1280 {{1, 2, 3, 4}}}};
1281
1282 return {&table[vertex][0], vertex == 4 ? 4u : 3u};
1283 }
1284 case ReferenceCells::Wedge:
1285 {
1286 AssertIndexRange(vertex, 6);
1287 static constexpr ndarray<unsigned int, 6, 3> table = {{{{0, 2, 4}},
1288 {{0, 2, 3}},
1289 {{0, 3, 4}},
1290 {{1, 2, 4}},
1291 {{1, 2, 3}},
1292 {{1, 3, 4}}}};
1293
1294 return table[vertex];
1295 }
1296 case ReferenceCells::Hexahedron:
1297 return {&GeometryInfo<3>::vertex_to_face[vertex][0], 3};
1298 default:
1299 DEAL_II_NOT_IMPLEMENTED();
1300 }
1301
1302 return {};
1303}
1304
1305
1306constexpr ::ndarray<unsigned int, 12, 4>
1307ReferenceCell::new_isotropic_child_face_lines(
1308 const unsigned int refienement_choice) const
1309{
1310 // AssertIndexRange(refienement_choice, n_isotropic_refinement_choices());
1311 const unsigned int X = numbers::invalid_unsigned_int;
1312 switch (this->kind)
1313 {
1314 case ReferenceCells::Tetrahedron:
1315 {
1316 constexpr ::ndarray<unsigned int, 3, 12, 4> new_quad_lines_tet =
1317 {{// new line is (6,8)
1318 {{{{2, 3, 8, X}},
1319 {{0, 9, 5, X}},
1320 {{1, 6, 11, X}},
1321 {{4, 10, 7, X}},
1322 {{2, 12, 5, X}},
1323 {{1, 9, 12, X}},
1324 {{4, 8, 12, X}},
1325 {{6, 12, 10, X}},
1326 {{X, X, X, X}},
1327 {{X, X, X, X}},
1328 {{X, X, X, X}},
1329 {{X, X, X, X}}}},
1330 // new line is (5,7)
1331 {{{{2, 3, 8, X}},
1332 {{0, 9, 5, X}},
1333 {{1, 6, 11, X}},
1334 {{4, 10, 7, X}},
1335 {{0, 3, 12, X}},
1336 {{1, 12, 8, X}},
1337 {{4, 12, 9, X}},
1338 {{7, 11, 12, X}},
1339 {{X, X, X, X}},
1340 {{X, X, X, X}},
1341 {{X, X, X, X}},
1342 {{X, X, X, X}}}},
1343 // new line is (4,9)
1344 {{{{2, 3, 8, X}},
1345 {{0, 9, 5, X}},
1346 {{1, 6, 11, X}},
1347 {{4, 10, 7, X}},
1348 {{0, 12, 11, X}},
1349 {{2, 6, 12, X}},
1350 {{3, 12, 7, X}},
1351 {{5, 10, 12, X}},
1352 {{X, X, X, X}},
1353 {{X, X, X, X}},
1354 {{X, X, X, X}},
1355 {{X, X, X, X}}}}}};
1356 return new_quad_lines_tet[refienement_choice];
1357 }
1358
1359 case ReferenceCells::Hexahedron:
1360 {
1361 constexpr ::ndarray<unsigned int, 12, 4> new_quad_lines_hex = {
1362 {{{10, 28, 16, 24}},
1363 {{28, 14, 17, 25}},
1364 {{11, 29, 24, 20}},
1365 {{29, 15, 25, 21}},
1366 {{18, 26, 0, 28}},
1367 {{26, 22, 1, 29}},
1368 {{19, 27, 28, 4}},
1369 {{27, 23, 29, 5}},
1370 {{2, 24, 8, 26}},
1371 {{24, 6, 9, 27}},
1372 {{3, 25, 26, 12}},
1373 {{25, 7, 27, 13}}}};
1374 return new_quad_lines_hex;
1375 }
1376 default:
1377 DEAL_II_NOT_IMPLEMENTED();
1378 }
1379
1380 return {};
1381}
1382
1383
1384
1385constexpr ::ndarray<unsigned int, 12, 4, 2>
1386ReferenceCell::new_isotropic_child_face_line_vertices(
1387 const unsigned int refienement_choice) const
1388{
1389 // AssertIndexRange(refienement_choice, n_isotropic_refinement_choices());
1390 const unsigned int X = numbers::invalid_unsigned_int;
1391 switch (this->kind)
1392 {
1393 case ReferenceCells::Tetrahedron:
1394 {
1395 constexpr ::ndarray<unsigned int, 3, 12, 4, 2> table_tet = {
1396 {// new line is (6, 8)
1397 {{{{{{6, 4}}, {{4, 7}}, {{7, 6}}, {{X, X}}}},
1398 {{{{4, 5}}, {{5, 8}}, {{8, 4}}, {{X, X}}}},
1399 {{{{5, 6}}, {{6, 9}}, {{9, 5}}, {{X, X}}}},
1400 {{{{7, 8}}, {{8, 9}}, {{9, 7}}, {{X, X}}}},
1401 {{{{4, 6}}, {{6, 8}}, {{8, 4}}, {{X, X}}}},
1402 {{{{6, 5}}, {{5, 8}}, {{8, 6}}, {{X, X}}}},
1403 {{{{8, 7}}, {{7, 6}}, {{6, 8}}, {{X, X}}}},
1404 {{{{9, 6}}, {{6, 8}}, {{8, 9}}, {{X, X}}}},
1405 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1406 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1407 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1408 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}}}},
1409 // new line is (5, 7)
1410 {{{{{{6, 4}}, {{4, 7}}, {{7, 6}}, {{X, X}}}},
1411 {{{{4, 5}}, {{5, 8}}, {{8, 4}}, {{X, X}}}},
1412 {{{{5, 6}}, {{6, 9}}, {{9, 5}}, {{X, X}}}},
1413 {{{{7, 8}}, {{8, 9}}, {{9, 7}}, {{X, X}}}},
1414 {{{{5, 4}}, {{4, 7}}, {{7, 5}}, {{X, X}}}},
1415 {{{{6, 5}}, {{5, 7}}, {{7, 6}}, {{X, X}}}},
1416 {{{{8, 7}}, {{7, 5}}, {{5, 8}}, {{X, X}}}},
1417 {{{{7, 9}}, {{9, 5}}, {{5, 7}}, {{X, X}}}},
1418 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1419 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1420 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1421 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}}}},
1422 // new line is (4, 9)
1423 {{{{{{6, 4}}, {{4, 7}}, {{7, 6}}, {{X, X}}}},
1424 {{{{4, 5}}, {{5, 8}}, {{8, 4}}, {{X, X}}}},
1425 {{{{5, 6}}, {{6, 9}}, {{9, 5}}, {{X, X}}}},
1426 {{{{7, 8}}, {{8, 9}}, {{9, 7}}, {{X, X}}}},
1427 {{{{5, 4}}, {{4, 9}}, {{9, 5}}, {{X, X}}}},
1428 {{{{4, 6}}, {{6, 9}}, {{9, 4}}, {{X, X}}}},
1429 {{{{7, 4}}, {{4, 9}}, {{9, 7}}, {{X, X}}}},
1430 {{{{4, 8}}, {{8, 9}}, {{9, 4}}, {{X, X}}}},
1431 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1432 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1433 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1434 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}}}}}};
1435 return table_tet[refienement_choice];
1436 }
1437
1438 case ReferenceCells::Hexahedron:
1439 {
1440 constexpr ::ndarray<unsigned int, 12, 4, 2> table_hex = {
1441 {{{{{10, 22}}, {{24, 26}}, {{10, 24}}, {{22, 26}}}},
1442 {{{{24, 26}}, {{11, 23}}, {{24, 11}}, {{26, 23}}}},
1443 {{{{22, 14}}, {{26, 25}}, {{22, 26}}, {{14, 25}}}},
1444 {{{{26, 25}}, {{23, 15}}, {{26, 23}}, {{25, 15}}}},
1445 {{{{8, 24}}, {{20, 26}}, {{8, 20}}, {{24, 26}}}},
1446 {{{{20, 26}}, {{12, 25}}, {{20, 12}}, {{26, 25}}}},
1447 {{{{24, 9}}, {{26, 21}}, {{24, 26}}, {{9, 21}}}},
1448 {{{{26, 21}}, {{25, 13}}, {{26, 25}}, {{21, 13}}}},
1449 {{{{16, 20}}, {{22, 26}}, {{16, 22}}, {{20, 26}}}},
1450 {{{{22, 26}}, {{17, 21}}, {{22, 17}}, {{26, 21}}}},
1451 {{{{20, 18}}, {{26, 23}}, {{20, 26}}, {{18, 23}}}},
1452 {{{{26, 23}}, {{21, 19}}, {{26, 21}}, {{23, 19}}}}}};
1453 return table_hex;
1454 }
1455 default:
1456 DEAL_II_NOT_IMPLEMENTED();
1457 }
1458
1459 return {};
1460}
1461
1462
1463
1464constexpr ::ndarray<unsigned int, 8, 6>
1465ReferenceCell::new_isotropic_child_cell_faces(
1466 const unsigned int refienement_choice) const
1467{
1468 // AssertIndexRange(refienement_choice, n_isotropic_refinement_choices());
1469 const unsigned int X = numbers::invalid_unsigned_int;
1470 switch (this->kind)
1471 {
1472 case ReferenceCells::Tetrahedron:
1473 {
1474 constexpr ::ndarray<unsigned int, 3, 8, 6> cell_quads_tet = {
1475 {// new line is (6, 8)
1476 {{{{8, 13, 16, 0, X, X}},
1477 {{9, 12, 1, 21, X, X}},
1478 {{10, 2, 17, 20, X, X}},
1479 {{3, 14, 18, 22, X, X}},
1480 {{11, 1, 4, 5, X, X}},
1481 {{15, 0, 4, 6, X, X}},
1482 {{19, 7, 6, 3, X, X}},
1483 {{23, 5, 2, 7, X, X}}}},
1484 // new line is (5, 7)
1485 {{
1486 {{8, 13, 16, 0, X, X}},
1487 {{9, 12, 1, 21, X, X}},
1488 {{10, 2, 17, 20, X, X}},
1489 {{3, 14, 18, 22, X, X}},
1490 {{11, 4, 0, 5, X, X}},
1491 {{15, 4, 1, 6, X, X}},
1492 {{19, 2, 5, 7, X, X}},
1493 {{23, 6, 7, 3, X, X}},
1494 }},
1495 // new line is (4, 9)
1496 {{
1497 {{8, 13, 16, 0, X, X}},
1498 {{9, 12, 1, 21, X, X}},
1499 {{10, 2, 17, 20, X, X}},
1500 {{3, 14, 18, 22, X, X}},
1501 {{11, 4, 5, 2, X, X}},
1502 {{15, 6, 7, 3, X, X}},
1503 {{19, 5, 0, 6, X, X}},
1504 {{23, 1, 4, 7, X, X}},
1505 }}}};
1506 return cell_quads_tet[refienement_choice];
1507 }
1508
1509 case ReferenceCells::Hexahedron:
1510 {
1511 constexpr ::ndarray<unsigned int, 8, 6> cell_quads_hex = {{
1512 {{12, 0, 20, 4, 28, 8}}, // bottom children
1513 {{0, 16, 22, 6, 29, 9}}, //
1514 {{13, 1, 4, 24, 30, 10}}, //
1515 {{1, 17, 6, 26, 31, 11}}, //
1516 {{14, 2, 21, 5, 8, 32}}, // top children
1517 {{2, 18, 23, 7, 9, 33}}, //
1518 {{15, 3, 5, 25, 10, 34}}, //
1519 {{3, 19, 7, 27, 11, 35}} //
1520 }};
1521 return cell_quads_hex;
1522 }
1523 default:
1524 DEAL_II_NOT_IMPLEMENTED();
1525 }
1526
1527 return {};
1528}
1529
1530
1531
1532constexpr ::ndarray<unsigned int, 8, 4>
1533ReferenceCell::new_isotropic_child_cell_vertices(
1534 const unsigned int refienement_choice) const
1535{
1536 // AssertIndexRange(refienement_choice, n_isotropic_refinement_choices());
1537
1538 switch (this->kind)
1539 {
1540 case ReferenceCells::Tetrahedron:
1541 {
1542 constexpr ::ndarray<unsigned int, 3, 8, 4> cell_vertices_tet = {
1543 {// new line is (6,8)
1544 {{{{0, 4, 6, 7}},
1545 {{4, 1, 5, 8}},
1546 {{6, 5, 2, 9}},
1547 {{7, 8, 9, 3}},
1548 {{4, 5, 6, 8}},
1549 {{4, 7, 8, 6}},
1550 {{6, 9, 7, 8}},
1551 {{5, 8, 9, 6}}}},
1552 // new line is (5,7)
1553 {{{{0, 4, 6, 7}},
1554 {{4, 1, 5, 8}},
1555 {{6, 5, 2, 9}},
1556 {{7, 8, 9, 3}},
1557 {{4, 5, 6, 7}},
1558 {{4, 7, 8, 5}},
1559 {{6, 9, 7, 5}},
1560 {{5, 8, 9, 7}}}},
1561 // new line is (4,9)
1562 {{{{0, 4, 6, 7}},
1563 {{4, 1, 5, 8}},
1564 {{6, 5, 2, 9}},
1565 {{7, 8, 9, 3}},
1566 {{4, 5, 6, 9}},
1567 {{4, 7, 8, 9}},
1568 {{6, 9, 7, 4}},
1569 {{5, 8, 9, 4}}}}}};
1570 return cell_vertices_tet[refienement_choice];
1571 }
1572 default:
1573 DEAL_II_NOT_IMPLEMENTED();
1574 }
1575
1576 return {};
1577}
1578
1579
1580
1581inline bool
1582ReferenceCell::is_hyper_cube() const
1583{
1584 return (*this == ReferenceCells::Vertex || *this == ReferenceCells::Line ||
1585 *this == ReferenceCells::Quadrilateral ||
1586 *this == ReferenceCells::Hexahedron);
1587}
1588
1589
1590
1591inline bool
1592ReferenceCell::is_simplex() const
1593{
1594 return (*this == ReferenceCells::Vertex || *this == ReferenceCells::Line ||
1595 *this == ReferenceCells::Triangle ||
1596 *this == ReferenceCells::Tetrahedron);
1597}
1598
1599
1600
1601inline unsigned int
1602ReferenceCell::get_dimension() const
1603{
1604 switch (this->kind)
1605 {
1606 case ReferenceCells::Vertex:
1607 return 0;
1608 case ReferenceCells::Line:
1609 return 1;
1610 case ReferenceCells::Triangle:
1611 case ReferenceCells::Quadrilateral:
1612 return 2;
1613 case ReferenceCells::Tetrahedron:
1614 case ReferenceCells::Pyramid:
1615 case ReferenceCells::Wedge:
1616 case ReferenceCells::Hexahedron:
1617 return 3;
1618 default:
1619 DEAL_II_NOT_IMPLEMENTED();
1620 }
1621
1622 return numbers::invalid_unsigned_int;
1623}
1624
1625
1626
1627template <int dim>
1628Quadrature<dim>
1629ReferenceCell::get_midpoint_quadrature() const
1630{
1631 Assert(dim == get_dimension(),
1632 ExcMessage("You can only call this function with arguments for "
1633 "which 'dim' equals the dimension of the space in which "
1634 "the current reference cell lives. You have dim=" +
1635 std::to_string(dim) + " but the reference cell is " +
1636 std::to_string(get_dimension()) + " dimensional."));
1637
1638 return Quadrature<dim>(std::vector<Point<dim>>({barycenter<dim>()}),
1639 std::vector<double>({volume()}));
1640}
1641
1642
1643
1644inline unsigned int
1645ReferenceCell::n_vertices() const
1646{
1647 switch (this->kind)
1648 {
1649 case ReferenceCells::Vertex:
1650 return 1;
1651 case ReferenceCells::Line:
1652 return 2;
1653 case ReferenceCells::Triangle:
1654 return 3;
1655 case ReferenceCells::Quadrilateral:
1656 return 4;
1657 case ReferenceCells::Tetrahedron:
1658 return 4;
1659 case ReferenceCells::Pyramid:
1660 return 5;
1661 case ReferenceCells::Wedge:
1662 return 6;
1663 case ReferenceCells::Hexahedron:
1664 return 8;
1665 default:
1666 DEAL_II_NOT_IMPLEMENTED();
1667 }
1668
1669 return numbers::invalid_unsigned_int;
1670}
1671
1672
1673
1674inline unsigned int
1675ReferenceCell::n_lines() const
1676{
1677 switch (this->kind)
1678 {
1679 case ReferenceCells::Vertex:
1680 return 0;
1681 case ReferenceCells::Line:
1682 return 1;
1683 case ReferenceCells::Triangle:
1684 return 3;
1685 case ReferenceCells::Quadrilateral:
1686 return 4;
1687 case ReferenceCells::Tetrahedron:
1688 return 6;
1689 case ReferenceCells::Pyramid:
1690 return 8;
1691 case ReferenceCells::Wedge:
1692 return 9;
1693 case ReferenceCells::Hexahedron:
1694 return 12;
1695 default:
1696 DEAL_II_NOT_IMPLEMENTED();
1697 }
1698
1699 return numbers::invalid_unsigned_int;
1700}
1701
1702
1703
1704template <int dim>
1705Point<dim>
1706ReferenceCell::vertex(const unsigned int v) const
1707{
1708 Assert(dim == get_dimension(),
1709 ExcMessage("You can only call this function with arguments for "
1710 "which 'dim' equals the dimension of the space in which "
1711 "the current reference cell lives. You have dim=" +
1712 std::to_string(dim) + " but the reference cell is " +
1713 std::to_string(get_dimension()) + " dimensional."));
1714 AssertIndexRange(v, n_vertices());
1715
1716 switch (dim)
1717 {
1718 case 0:
1719 {
1720 if (*this == ReferenceCells::Vertex)
1721 return Point<dim>();
1722 break;
1723 }
1724 case 1:
1725 {
1726 static const Point<dim> vertices[2] = {
1727 Point<dim>(), // the origin
1728 Point<dim>::unit_vector(0) // unit point along x-axis
1729 };
1730 if (*this == ReferenceCells::Line)
1731 return vertices[v];
1732 break;
1733 }
1734 case 2:
1735 {
1736 switch (this->kind)
1737 {
1738 case ReferenceCells::Triangle:
1739 {
1740 static const Point<dim> vertices[3] = {
1741 Point<dim>(), // the origin
1742 Point<dim>::unit_vector(0), // unit point along x-axis
1743 Point<dim>::unit_vector(1) // unit point along y-axis
1744 };
1745 return vertices[v];
1746 }
1747 case ReferenceCells::Quadrilateral:
1748 {
1749 static const Point<dim> vertices[4] = {
1750 // First the two points on the x-axis
1751 Point<dim>(),
1752 Point<dim>::unit_vector(0),
1753 // Then these two points shifted in the y-direction
1754 Point<dim>() + Point<dim>::unit_vector(1),
1755 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(1)};
1756 return vertices[v];
1757 }
1758 }
1759 break;
1760 }
1761 case 3:
1762 {
1763 switch (this->kind)
1764 {
1765 case ReferenceCells::Tetrahedron:
1766 {
1767 static const Point<dim> vertices[4] = {
1768 Point<dim>(), // the origin
1769 Point<dim>::unit_vector(0), // unit point along x-axis
1770 Point<dim>::unit_vector(1), // unit point along y-axis
1771 Point<dim>::unit_vector(2) // unit point along z-axis
1772 };
1773 return vertices[v];
1774 }
1775 case ReferenceCells::Pyramid:
1776 {
1777 static const Point<dim> vertices[5] = {
1778 Point<dim>{-1.0, -1.0, 0.0},
1779 Point<dim>{+1.0, -1.0, 0.0},
1780 Point<dim>{-1.0, +1.0, 0.0},
1781 Point<dim>{+1.0, +1.0, 0.0},
1782 Point<dim>{+0.0, +0.0, 1.0}};
1783 return vertices[v];
1784 }
1785 case ReferenceCells::Wedge:
1786 {
1787 static const Point<dim> vertices[6] = {
1788 // First the three points on the triangular base of the
1789 // wedge:
1790 Point<dim>(),
1791 Point<dim>::unit_vector(0),
1792 Point<dim>::unit_vector(1),
1793 // And now everything shifted in the z-direction again
1794 Point<dim>() + Point<dim>::unit_vector(2),
1795 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(2),
1796 Point<dim>::unit_vector(1) + Point<dim>::unit_vector(2)};
1797 return vertices[v];
1798 }
1799 case ReferenceCells::Hexahedron:
1800 {
1801 static const Point<dim> vertices[8] = {
1802 // First the two points on the x-axis
1803 Point<dim>(),
1804 Point<dim>::unit_vector(0),
1805 // Then these two points shifted in the y-direction
1806 Point<dim>() + Point<dim>::unit_vector(1),
1807 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(1),
1808 // And now all four points shifted in the z-direction
1809 Point<dim>() + Point<dim>::unit_vector(2),
1810 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(2),
1811 Point<dim>() + Point<dim>::unit_vector(1) +
1812 Point<dim>::unit_vector(2),
1813 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(1) +
1814 Point<dim>::unit_vector(2)};
1815 return vertices[v];
1816 }
1817 }
1818 break;
1819 }
1820 default:
1821 DEAL_II_NOT_IMPLEMENTED();
1822 }
1823
1824 DEAL_II_NOT_IMPLEMENTED();
1825 return Point<dim>();
1826}
1827
1828
1829inline unsigned int
1830ReferenceCell::n_faces() const
1831{
1832 switch (this->kind)
1833 {
1834 case ReferenceCells::Vertex:
1835 return 0;
1836 case ReferenceCells::Line:
1837 return 2;
1838 case ReferenceCells::Triangle:
1839 return 3;
1840 case ReferenceCells::Quadrilateral:
1841 return 4;
1842 case ReferenceCells::Tetrahedron:
1843 return 4;
1844 case ReferenceCells::Pyramid:
1845 return 5;
1846 case ReferenceCells::Wedge:
1847 return 5;
1848 case ReferenceCells::Hexahedron:
1849 return 6;
1850 default:
1851 DEAL_II_NOT_IMPLEMENTED();
1852 }
1853
1854 return numbers::invalid_unsigned_int;
1855}
1856
1857
1858
1859inline std_cxx20::ranges::iota_view<unsigned int, unsigned int>
1860ReferenceCell::face_indices() const
1861{
1862 return {0U, n_faces()};
1863}
1864
1865
1866
1867inline unsigned int
1868ReferenceCell::n_isotropic_refinement_choices() const
1869{
1870 switch (this->kind)
1871 {
1872 case ReferenceCells::Vertex:
1873 return 1;
1874 case ReferenceCells::Line:
1875 return 1;
1876 case ReferenceCells::Triangle:
1877 return 1;
1878 case ReferenceCells::Quadrilateral:
1879 return 1;
1880 case ReferenceCells::Tetrahedron:
1881 return 3;
1882 case ReferenceCells::Pyramid:
1883 return 1;
1884 case ReferenceCells::Wedge:
1885 return 1;
1886 case ReferenceCells::Hexahedron:
1887 return 1;
1888 default:
1889 DEAL_II_ASSERT_UNREACHABLE();
1890 }
1891
1892 return numbers::invalid_unsigned_int;
1893}
1894
1895
1896
1897inline IsotropicRefinementChoice
1898ReferenceCell::get_isotropic_refinement_choice(
1899 const unsigned int ref_choice) const
1900{
1901 AssertIndexRange(ref_choice, n_isotropic_refinement_choices());
1902 switch (this->kind)
1903 {
1904 case ReferenceCells::Tetrahedron:
1905 {
1906 static constexpr ndarray<IsotropicRefinementChoice, 3>
1907 isotropic_ref_choices = {{IsotropicRefinementChoice::cut_tet_68,
1908 IsotropicRefinementChoice::cut_tet_57,
1909 IsotropicRefinementChoice::cut_tet_49}};
1910 return isotropic_ref_choices[ref_choice];
1911 }
1912 default:
1913 {
1914 DEAL_II_NOT_IMPLEMENTED();
1915 }
1916 }
1917
1918 return IsotropicRefinementChoice::isotropic_refinement;
1919}
1920
1921
1922
1923inline unsigned int
1924ReferenceCell::n_isotropic_children() const
1925{
1926 switch (this->kind)
1927 {
1928 case ReferenceCells::Vertex:
1929 return 0;
1930 case ReferenceCells::Line:
1931 return 2;
1932 case ReferenceCells::Triangle:
1933 return 4;
1934 case ReferenceCells::Quadrilateral:
1935 return 4;
1936 case ReferenceCells::Tetrahedron:
1937 return 8;
1938 case ReferenceCells::Pyramid:
1939 // We haven't yet decided how to refine pyramids. Update this when we
1940 // have
1941 return 0;
1942 case ReferenceCells::Wedge:
1943 return 8;
1944 case ReferenceCells::Hexahedron:
1945 return 8;
1946 default:
1947 DEAL_II_ASSERT_UNREACHABLE();
1948 }
1949
1950 return numbers::invalid_unsigned_int;
1951}
1952
1953
1954
1955template <int dim>
1956unsigned int
1957ReferenceCell::n_children(const RefinementCase<dim> ref_case) const
1958{
1959 Assert(dim == get_dimension(),
1960 ExcMessage("You can only call this function with arguments for "
1961 "which 'dim' equals the dimension of the space in which "
1962 "the current reference cell lives. You have dim=" +
1963 std::to_string(dim) + " but the reference cell is " +
1964 std::to_string(get_dimension()) + " dimensional."));
1965
1966 // Use GeometryInfo here to keep it the single source of truth
1967 if (this->is_hyper_cube())
1968 return GeometryInfo<dim>::n_children(ref_case);
1969 else
1970 return this->n_isotropic_children();
1971}
1972
1973
1974
1975inline std_cxx20::ranges::iota_view<unsigned int, unsigned int>
1976ReferenceCell::isotropic_child_indices() const
1977{
1978 return {0U, n_isotropic_children()};
1979}
1980
1981
1982
1983inline std_cxx20::ranges::iota_view<unsigned int, unsigned int>
1984ReferenceCell::vertex_indices() const
1985{
1986 return {0U, n_vertices()};
1987}
1988
1989
1990
1991inline std_cxx20::ranges::iota_view<unsigned int, unsigned int>
1992ReferenceCell::line_indices() const
1993{
1994 return {0U, n_lines()};
1995}
1996
1997
1998
1999inline ReferenceCell
2000ReferenceCell::face_reference_cell(const unsigned int face_no) const
2001{
2002 AssertIndexRange(face_no, n_faces());
2003
2004 switch (this->kind)
2005 {
2006 case ReferenceCells::Vertex:
2007 return ReferenceCells::Invalid;
2008 case ReferenceCells::Line:
2009 return ReferenceCells::Vertex;
2010 case ReferenceCells::Triangle:
2011 case ReferenceCells::Quadrilateral:
2012 return ReferenceCells::Line;
2013 case ReferenceCells::Tetrahedron:
2014 return ReferenceCells::Triangle;
2015 case ReferenceCells::Pyramid:
2016 if (face_no == 0)
2017 return ReferenceCells::Quadrilateral;
2018 else
2019 return ReferenceCells::Triangle;
2020 case ReferenceCells::Wedge:
2021 if (face_no > 1)
2022 return ReferenceCells::Quadrilateral;
2023 else
2024 return ReferenceCells::Triangle;
2025 case ReferenceCells::Hexahedron:
2026 return ReferenceCells::Quadrilateral;
2027 default:
2028 DEAL_II_NOT_IMPLEMENTED();
2029 }
2030
2031 return ReferenceCells::Invalid;
2032}
2033
2034
2035
2036inline constexpr types::geometric_orientation
2037ReferenceCell::default_combined_face_orientation()
2038{
2039 // Our convention is that 'orientation' has a default value of true and
2040 // occupies the least-significant bit while rotate and flip have default
2041 // values of 'false' and occupy the second and third bits.
2042 return 0b001;
2043}
2044
2045
2046
2047inline constexpr types::geometric_orientation
2048ReferenceCell::reversed_combined_line_orientation()
2049{
2050 // For a reversed line 'orientation' is false and neither flip nor rotate are
2051 // defined.
2052 return 0b000;
2053}
2054
2055
2056
2057inline unsigned int
2058ReferenceCell::child_cell_on_face(
2059 const unsigned int face,
2060 const unsigned int subface,
2061 const types::geometric_orientation combined_face_orientation) const
2062{
2063 AssertIndexRange(face, n_faces());
2064 AssertIndexRange(subface, face_reference_cell(face).n_isotropic_children());
2065
2066 switch (this->kind)
2067 {
2068 case ReferenceCells::Vertex:
2069 case ReferenceCells::Line:
2070 {
2071 DEAL_II_NOT_IMPLEMENTED();
2072 break;
2073 }
2074 case ReferenceCells::Triangle:
2075 {
2076 static constexpr ndarray<unsigned int, 3, 2> subcells = {
2077 {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
2078
2079 return subcells[face][subface];
2080 }
2081 case ReferenceCells::Quadrilateral:
2082 {
2083 const auto [face_orientation, face_rotation, face_flip] =
2084 internal::split_face_orientation(combined_face_orientation);
2085
2086 return GeometryInfo<2>::child_cell_on_face(
2087 RefinementCase<2>(RefinementPossibilities<2>::isotropic_refinement),
2088 face,
2089 subface,
2090 face_orientation,
2091 face_flip,
2092 face_rotation);
2093 }
2094 case ReferenceCells::Tetrahedron:
2095 case ReferenceCells::Pyramid:
2096 case ReferenceCells::Wedge:
2097 {
2098 DEAL_II_NOT_IMPLEMENTED();
2099 break;
2100 }
2101 case ReferenceCells::Hexahedron:
2102 {
2103 const auto [face_orientation, face_rotation, face_flip] =
2104 internal::split_face_orientation(combined_face_orientation);
2105
2106 return GeometryInfo<3>::child_cell_on_face(
2107 RefinementCase<3>(RefinementPossibilities<3>::isotropic_refinement),
2108 face,
2109 subface,
2110 face_orientation,
2111 face_flip,
2112 face_rotation);
2113 }
2114 default:
2115 DEAL_II_NOT_IMPLEMENTED();
2116 }
2117
2118 return numbers::invalid_unsigned_int;
2119}
2120
2121
2122
2123inline std::array<unsigned int, 2>
2124ReferenceCell::standard_vertex_to_face_and_vertex_index(
2125 const unsigned int vertex) const
2126{
2127 AssertIndexRange(vertex, n_vertices());
2128 // Work around a GCC warning at higher optimization levels by making all of
2129 // these tables the same size
2130 constexpr unsigned int X = numbers::invalid_unsigned_int;
2131
2132 switch (this->kind)
2133 {
2134 case ReferenceCells::Vertex:
2135 case ReferenceCells::Line:
2136 DEAL_II_NOT_IMPLEMENTED();
2137 break;
2138 case ReferenceCells::Triangle:
2139 {
2140 static constexpr ndarray<unsigned int, 6, 2> table = {
2141 {{{0, 0}}, {{0, 1}}, {{1, 1}}, {{X, X}}, {{X, X}}, {{X, X}}}};
2142
2143 return table[vertex];
2144 }
2145 case ReferenceCells::Quadrilateral:
2146 {
2147 return GeometryInfo<2>::standard_quad_vertex_to_line_vertex_index(
2148 vertex);
2149 }
2150 case ReferenceCells::Tetrahedron:
2151 {
2152 static constexpr ndarray<unsigned int, 6, 2> table = {
2153 {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 2}}, {{X, X}}, {{X, X}}}};
2154
2155 return table[vertex];
2156 }
2157 case ReferenceCells::Pyramid:
2158 {
2159 static constexpr ndarray<unsigned int, 6, 2> table = {
2160 {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{0, 3}}, {{1, 2}}, {{X, X}}}};
2161
2162 return table[vertex];
2163 }
2164 case ReferenceCells::Wedge:
2165 {
2166 static constexpr ndarray<unsigned int, 6, 2> table = {
2167 {{{0, 1}}, {{0, 0}}, {{0, 2}}, {{1, 0}}, {{1, 1}}, {{1, 2}}}};
2168
2169 return table[vertex];
2170 }
2171 case ReferenceCells::Hexahedron:
2172 {
2173 return GeometryInfo<3>::standard_hex_vertex_to_quad_vertex_index(
2174 vertex);
2175 }
2176 default:
2177 DEAL_II_NOT_IMPLEMENTED();
2178 }
2179
2180 return {};
2181}
2182
2183
2184
2185inline std::array<unsigned int, 2>
2186ReferenceCell::standard_line_to_face_and_line_index(
2187 const unsigned int line) const
2188{
2189 AssertIndexRange(line, n_lines());
2190
2191 switch (this->kind)
2192 {
2193 case ReferenceCells::Vertex:
2194 case ReferenceCells::Line:
2195 case ReferenceCells::Triangle:
2196 case ReferenceCells::Quadrilateral:
2197 {
2198 DEAL_II_NOT_IMPLEMENTED();
2199 break;
2200 }
2201 case ReferenceCells::Tetrahedron:
2202 {
2203 static const std::array<unsigned int, 2> table[6] = {
2204 {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}};
2205
2206 return table[line];
2207 }
2208 case ReferenceCells::Pyramid:
2209 {
2210 static const std::array<unsigned int, 2> table[8] = {{{0, 0}},
2211 {{0, 1}},
2212 {{0, 2}},
2213 {{0, 3}},
2214 {{1, 2}},
2215 {{2, 1}},
2216 {{1, 1}},
2217 {{2, 2}}};
2218
2219 return table[line];
2220 }
2221 case ReferenceCells::Wedge:
2222 {
2223 static const std::array<unsigned int, 2> table[9] = {{{0, 0}},
2224 {{0, 2}},
2225 {{0, 1}},
2226 {{1, 0}},
2227 {{1, 1}},
2228 {{1, 2}},
2229 {{2, 0}},
2230 {{2, 1}},
2231 {{3, 1}}};
2232
2233 return table[line];
2234 }
2235 case ReferenceCells::Hexahedron:
2236 {
2237 return GeometryInfo<3>::standard_hex_line_to_quad_line_index(line);
2238 }
2239 default:
2240 DEAL_II_NOT_IMPLEMENTED();
2241 }
2242
2243 return {};
2244}
2245
2246
2247inline unsigned int
2248ReferenceCell::line_to_cell_vertices(const unsigned int line,
2249 const unsigned int vertex) const
2250{
2251 AssertIndexRange(vertex, 2);
2252 AssertIndexRange(line, n_lines());
2253
2254 switch (this->kind)
2255 {
2256 case ReferenceCells::Vertex:
2257 case ReferenceCells::Line:
2258 return vertex;
2259 case ReferenceCells::Triangle:
2260 {
2261 static constexpr ndarray<unsigned int, 3, 2> table = {
2262 {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
2263 return table[line][vertex];
2264 }
2265 case ReferenceCells::Quadrilateral:
2266 {
2267 static constexpr ndarray<unsigned int, 4, 2> table = {
2268 {{{0, 2}}, {{1, 3}}, {{0, 1}}, {{2, 3}}}};
2269 return table[line][vertex];
2270 }
2271 case ReferenceCells::Tetrahedron:
2272 {
2273 static constexpr ndarray<unsigned int, 6, 2> table = {
2274 {{{0, 1}}, {{1, 2}}, {{2, 0}}, {{0, 3}}, {{1, 3}}, {{2, 3}}}};
2275 return table[line][vertex];
2276 }
2277 case ReferenceCells::Pyramid:
2278 {
2279 static constexpr ndarray<unsigned int, 8, 2> table = {{{{0, 2}},
2280 {{1, 3}},
2281 {{0, 1}},
2282 {{2, 3}},
2283 {{4, 0}},
2284 {{1, 4}},
2285 {{2, 4}},
2286 {{4, 3}}}};
2287 return table[line][vertex];
2288 }
2289 case ReferenceCells::Wedge:
2290 {
2291 static constexpr ndarray<unsigned int, 9, 2> table = {{{{1, 0}},
2292 {{2, 1}},
2293 {{0, 2}},
2294 {{3, 4}},
2295 {{4, 5}},
2296 {{5, 3}},
2297 {{0, 3}},
2298 {{1, 4}},
2299 {{2, 5}}}};
2300 return table[line][vertex];
2301 }
2302 case ReferenceCells::Hexahedron:
2303 {
2304 // first four lines comprise the bottom face, next four are the top,
2305 // and the last four are 'bottom to top'
2306 static constexpr ndarray<unsigned int, 12, 2> table = {{{{0, 2}},
2307 {{1, 3}},
2308 {{0, 1}},
2309 {{2, 3}},
2310 {{4, 6}},
2311 {{5, 7}},
2312 {{4, 5}},
2313 {{6, 7}},
2314 {{0, 4}},
2315 {{1, 5}},
2316 {{2, 6}},
2317 {{3, 7}}}};
2318 return table[line][vertex];
2319 }
2320
2321 default:
2322 DEAL_II_NOT_IMPLEMENTED();
2323 }
2324
2325 return numbers::invalid_unsigned_int;
2326}
2327
2328
2329inline unsigned int
2330ReferenceCell::face_to_cell_lines(
2331 const unsigned int face,
2332 const unsigned int line,
2333 const types::geometric_orientation combined_face_orientation) const
2334{
2335 AssertIndexRange(face, n_faces());
2336 AssertIndexRange(line, face_reference_cell(face).n_lines());
2337
2338 static constexpr unsigned int X = numbers::invalid_unsigned_int;
2339
2340 switch (this->kind)
2341 {
2342 case ReferenceCells::Vertex:
2343 {
2344 DEAL_II_NOT_IMPLEMENTED();
2345 break;
2346 }
2347 case ReferenceCells::Line:
2348 {
2349 const auto [face_orientation, face_rotation, face_flip] =
2350 internal::split_face_orientation(combined_face_orientation);
2351
2352 return GeometryInfo<1>::face_to_cell_lines(
2353 face, line, face_orientation, face_flip, face_rotation);
2354 }
2355 case ReferenceCells::Triangle:
2356 {
2357 return face;
2358 }
2359 case ReferenceCells::Quadrilateral:
2360 {
2361 const auto [face_orientation, face_rotation, face_flip] =
2362 internal::split_face_orientation(combined_face_orientation);
2363
2364 return GeometryInfo<2>::face_to_cell_lines(
2365 face, line, face_orientation, face_flip, face_rotation);
2366 }
2367 case ReferenceCells::Tetrahedron:
2368 {
2369 static constexpr ndarray<unsigned int, 4, 3> table = {
2370 {{{0, 1, 2}}, {{0, 3, 4}}, {{2, 5, 3}}, {{1, 4, 5}}}};
2371
2372 return table[face][standard_to_real_face_line(
2373 line, face, combined_face_orientation)];
2374 }
2375 case ReferenceCells::Pyramid:
2376 {
2377 static constexpr ndarray<unsigned int, 5, 4> table = {
2378 {{{0, 1, 2, 3}},
2379 {{0, 6, 4, X}},
2380 {{1, 5, 7, X}},
2381 {{2, 4, 5, X}},
2382 {{3, 7, 6, X}}}};
2383
2384 return table[face][standard_to_real_face_line(
2385 line, face, combined_face_orientation)];
2386 }
2387 case ReferenceCells::Wedge:
2388 {
2389 static constexpr ndarray<unsigned int, 5, 4> table = {
2390 {{{0, 2, 1, X}},
2391 {{3, 4, 5, X}},
2392 {{6, 7, 0, 3}},
2393 {{7, 8, 1, 4}},
2394 {{8, 6, 5, 2}}}};
2395
2396 return table[face][standard_to_real_face_line(
2397 line, face, combined_face_orientation)];
2398 }
2399 case ReferenceCells::Hexahedron:
2400 {
2401 const auto [face_orientation, face_rotation, face_flip] =
2402 internal::split_face_orientation(combined_face_orientation);
2403
2404 return GeometryInfo<3>::face_to_cell_lines(
2405 face, line, face_orientation, face_flip, face_rotation);
2406 }
2407 default:
2408 DEAL_II_NOT_IMPLEMENTED();
2409 }
2410
2411 return numbers::invalid_unsigned_int;
2412}
2413
2414
2415
2416inline unsigned int
2417ReferenceCell::face_to_cell_vertices(
2418 const unsigned int face,
2419 const unsigned int vertex,
2420 const types::geometric_orientation combined_face_orientation) const
2421{
2422 AssertIndexRange(face, n_faces());
2423 AssertIndexRange(vertex, face_reference_cell(face).n_vertices());
2424 // TODO: once the default orientation is switched to 0 then we can remove this
2425 // special case for 1D.
2426 if (get_dimension() == 1)
2427 Assert(combined_face_orientation ==
2428 ReferenceCell::default_combined_face_orientation(),
2429 ExcMessage("In 1D, all faces must have the default orientation."));
2430 else
2431 AssertIndexRange(combined_face_orientation, n_face_orientations(face));
2432
2433 switch (this->kind)
2434 {
2435 case ReferenceCells::Vertex:
2436 {
2437 DEAL_II_NOT_IMPLEMENTED();
2438 break;
2439 }
2440 case ReferenceCells::Line:
2441 {
2442 const auto [face_orientation, face_rotation, face_flip] =
2443 internal::split_face_orientation(combined_face_orientation);
2444
2445 return GeometryInfo<1>::face_to_cell_vertices(
2446 face, vertex, face_orientation, face_flip, face_rotation);
2447 }
2448 case ReferenceCells::Triangle:
2449 {
2450 static constexpr ndarray<unsigned int, 3, 2> table = {
2451 {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
2452
2453 return table[face][combined_face_orientation ==
2454 default_combined_face_orientation() ?
2455 vertex :
2456 (1 - vertex)];
2457 }
2458 case ReferenceCells::Quadrilateral:
2459 {
2460 const auto [face_orientation, face_rotation, face_flip] =
2461 internal::split_face_orientation(combined_face_orientation);
2462
2463 return GeometryInfo<2>::face_to_cell_vertices(
2464 face, vertex, face_orientation, face_flip, face_rotation);
2465 }
2466 case ReferenceCells::Tetrahedron:
2467 {
2468 static constexpr ndarray<unsigned int, 4, 3> table = {
2469 {{{0, 1, 2}}, {{1, 0, 3}}, {{0, 2, 3}}, {{2, 1, 3}}}};
2470
2471 return table[face][standard_to_real_face_vertex(
2472 vertex, face, combined_face_orientation)];
2473 }
2474 case ReferenceCells::Pyramid:
2475 {
2476 constexpr auto X = numbers::invalid_unsigned_int;
2477 static constexpr ndarray<unsigned int, 5, 4> table = {
2478 {{{0, 1, 2, 3}},
2479 {{0, 2, 4, X}},
2480 {{3, 1, 4, X}},
2481 {{1, 0, 4, X}},
2482 {{2, 3, 4, X}}}};
2483
2484 return table[face][standard_to_real_face_vertex(
2485 vertex, face, combined_face_orientation)];
2486 }
2487 case ReferenceCells::Wedge:
2488 {
2489 constexpr auto X = numbers::invalid_unsigned_int;
2490 static constexpr ndarray<unsigned int, 6, 4> table = {
2491 {{{1, 0, 2, X}},
2492 {{3, 4, 5, X}},
2493 {{0, 1, 3, 4}},
2494 {{1, 2, 4, 5}},
2495 {{2, 0, 5, 3}}}};
2496
2497 return table[face][standard_to_real_face_vertex(
2498 vertex, face, combined_face_orientation)];
2499 }
2500 case ReferenceCells::Hexahedron:
2501 {
2502 const auto [face_orientation, face_rotation, face_flip] =
2503 internal::split_face_orientation(combined_face_orientation);
2504
2505 return GeometryInfo<3>::face_to_cell_vertices(
2506 face, vertex, face_orientation, face_flip, face_rotation);
2507 }
2508 default:
2509 DEAL_II_NOT_IMPLEMENTED();
2510 }
2511
2512 return numbers::invalid_unsigned_int;
2513}
2514
2515
2516
2517template <int dim>
2518Point<dim>
2519ReferenceCell::face_vertex_location(const unsigned int face,
2520 const unsigned int vertex) const
2521{
2522 Assert(dim == get_dimension(),
2523 ExcMessage("You can only call this function with arguments for "
2524 "which 'dim' equals the dimension of the space in which "
2525 "the current reference cell lives. You have dim=" +
2526 std::to_string(dim) + " but the reference cell is " +
2527 std::to_string(get_dimension()) + " dimensional."));
2528 return this->template vertex<dim>(
2529 face_to_cell_vertices(face, vertex, default_combined_face_orientation()));
2530}
2531
2532
2533
2534inline unsigned int
2535ReferenceCell::standard_to_real_face_vertex(
2536 const unsigned int vertex,
2537 const unsigned int face,
2538 const types::geometric_orientation face_orientation) const
2539{
2540 AssertIndexRange(face, n_faces());
2541 AssertIndexRange(vertex, face_reference_cell(face).n_vertices());
2542
2543 switch (this->kind)
2544 {
2545 case ReferenceCells::Vertex:
2546 DEAL_II_NOT_IMPLEMENTED();
2547 break;
2548 case ReferenceCells::Line:
2549 Assert(face_orientation == default_combined_face_orientation(),
2550 ExcMessage(
2551 "In 1D, all faces must have the default orientation."));
2552 return vertex;
2553 case ReferenceCells::Triangle:
2554 case ReferenceCells::Quadrilateral:
2555 return line_vertex_permutations[face_orientation][vertex];
2556 case ReferenceCells::Tetrahedron:
2557 return triangle_vertex_permutations[face_orientation][vertex];
2558 case ReferenceCells::Pyramid:
2559 // face 0 is a quadrilateral
2560 if (face == 0)
2561 return quadrilateral_vertex_permutations[face_orientation][vertex];
2562 else
2563 return triangle_vertex_permutations[face_orientation][vertex];
2564 case ReferenceCells::Wedge:
2565 // faces 0 and 1 are triangles
2566 if (face > 1)
2567 return quadrilateral_vertex_permutations[face_orientation][vertex];
2568 else
2569 return triangle_vertex_permutations[face_orientation][vertex];
2570 case ReferenceCells::Hexahedron:
2571 return quadrilateral_vertex_permutations[face_orientation][vertex];
2572 default:
2573 DEAL_II_NOT_IMPLEMENTED();
2574 }
2575
2576 DEAL_II_NOT_IMPLEMENTED();
2577 return numbers::invalid_unsigned_int;
2578}
2579
2580
2581
2582inline unsigned int
2583ReferenceCell::standard_to_real_face_line(
2584 const unsigned int line,
2585 const unsigned int face,
2586 const types::geometric_orientation combined_face_orientation) const
2587{
2588 AssertIndexRange(face, n_faces());
2589 AssertIndexRange(line, face_reference_cell(face).n_lines());
2590
2591 static constexpr ndarray<unsigned int, 6, 3> triangle_table = {{{{2, 1, 0}},
2592 {{0, 1, 2}},
2593 {{1, 0, 2}},
2594 {{2, 0, 1}},
2595 {{0, 2, 1}},
2596 {{1, 2, 0}}}};
2597
2598
2599 switch (this->kind)
2600 {
2601 case ReferenceCells::Vertex:
2602 case ReferenceCells::Line:
2603 case ReferenceCells::Triangle:
2604 case ReferenceCells::Quadrilateral:
2605 DEAL_II_NOT_IMPLEMENTED();
2606 break;
2607 case ReferenceCells::Tetrahedron:
2608 return triangle_table[combined_face_orientation][line];
2609 case ReferenceCells::Pyramid:
2610 if (face == 0) // The quadrilateral face
2611 {
2612 const auto [face_orientation, face_rotation, face_flip] =
2613 internal::split_face_orientation(combined_face_orientation);
2614
2615 return GeometryInfo<3>::standard_to_real_face_line(line,
2616 face_orientation,
2617 face_flip,
2618 face_rotation);
2619 }
2620 else // One of the triangular faces
2621 {
2622 return triangle_table[combined_face_orientation][line];
2623 }
2624 case ReferenceCells::Wedge:
2625 if (face > 1) // One of the quadrilateral faces
2626 {
2627 const auto [face_orientation, face_rotation, face_flip] =
2628 internal::split_face_orientation(combined_face_orientation);
2629
2630 return GeometryInfo<3>::standard_to_real_face_line(line,
2631 face_orientation,
2632 face_flip,
2633 face_rotation);
2634 }
2635 else // One of the triangular faces
2636 return triangle_table[combined_face_orientation][line];
2637 case ReferenceCells::Hexahedron:
2638 {
2639 static constexpr ndarray<unsigned int, 8, 4> table = {
2640 {{{2, 3, 0, 1}},
2641 {{0, 1, 2, 3}},
2642 {{0, 1, 3, 2}},
2643 {{3, 2, 0, 1}},
2644 {{3, 2, 1, 0}},
2645 {{1, 0, 3, 2}},
2646 {{1, 0, 2, 3}},
2647 {{2, 3, 1, 0}}}};
2648 return table[combined_face_orientation][line];
2649 }
2650 default:
2651 DEAL_II_NOT_IMPLEMENTED();
2652 }
2653
2654 return numbers::invalid_unsigned_int;
2655}
2656
2657
2658
2659namespace ReferenceCells
2660{
2661 template <int dim>
2662 inline constexpr const ReferenceCell &
2663 get_simplex()
2664 {
2665 switch (dim)
2666 {
2667 case 0:
2668 return ReferenceCells::Vertex;
2669 case 1:
2670 return ReferenceCells::Line;
2671 case 2:
2672 return ReferenceCells::Triangle;
2673 case 3:
2674 return ReferenceCells::Tetrahedron;
2675 default:
2676 DEAL_II_NOT_IMPLEMENTED();
2677 }
2678 return ReferenceCells::Invalid;
2679 }
2680
2681
2682
2683 template <int dim>
2684 inline constexpr const ReferenceCell &
2685 get_hypercube()
2686 {
2687 switch (dim)
2688 {
2689 case 0:
2690 return ReferenceCells::Vertex;
2691 case 1:
2692 return ReferenceCells::Line;
2693 case 2:
2694 return ReferenceCells::Quadrilateral;
2695 case 3:
2696 return ReferenceCells::Hexahedron;
2697 default:
2698 DEAL_II_NOT_IMPLEMENTED();
2699 }
2700 return ReferenceCells::Invalid;
2701 }
2702
2703
2704
2705 template <int structdim>
2706 inline constexpr unsigned int
2711
2712
2713
2714 template <int structdim>
2715 inline constexpr unsigned int
2720
2721
2722
2723 template <int structdim>
2724 inline constexpr unsigned int
2729} // namespace ReferenceCells
2730
2731
2732inline ReferenceCell
2733ReferenceCell::n_vertices_to_type(const int dim, const unsigned int n_vertices)
2734{
2735 AssertIndexRange(dim, 4);
2736 AssertIndexRange(n_vertices, 9);
2737
2738 const auto X = ReferenceCells::Invalid;
2739 static const ndarray<ReferenceCell, 4, 9> table = {
2740 {// dim 0
2741 {{X, ReferenceCells::Vertex, X, X, X, X, X, X, X}},
2742 // dim 1
2743 {{X, X, ReferenceCells::Line, X, X, X, X, X, X}},
2744 // dim 2
2745 {{X,
2746 X,
2747 X,
2748 ReferenceCells::Triangle,
2749 ReferenceCells::Quadrilateral,
2750 X,
2751 X,
2752 X,
2753 X}},
2754 // dim 3
2755 {{X,
2756 X,
2757 X,
2758 X,
2759 ReferenceCells::Tetrahedron,
2760 ReferenceCells::Pyramid,
2761 ReferenceCells::Wedge,
2762 X,
2763 ReferenceCells::Hexahedron}}}};
2764 Assert(table[dim][n_vertices] != ReferenceCells::Invalid,
2765 ExcMessage("The combination of dim = " + std::to_string(dim) +
2766 " and n_vertices = " + std::to_string(n_vertices) +
2767 " does not correspond to a known reference cell type."));
2768 return table[dim][n_vertices];
2769}
2770
2771
2772
2773template <int dim>
2774inline double
2775ReferenceCell::d_linear_shape_function(const Point<dim> &xi,
2776 const unsigned int i) const
2777{
2778 Assert(dim == get_dimension(),
2779 ExcMessage("You can only call this function with arguments for "
2780 "which 'dim' equals the dimension of the space in which "
2781 "the current reference cell lives. You have dim=" +
2782 std::to_string(dim) + " but the reference cell is " +
2783 std::to_string(get_dimension()) + " dimensional."));
2784 AssertIndexRange(i, n_vertices());
2785 switch (this->kind)
2786 {
2787 case ReferenceCells::Vertex:
2788 case ReferenceCells::Line:
2789 case ReferenceCells::Quadrilateral:
2790 case ReferenceCells::Hexahedron:
2791 return GeometryInfo<dim>::d_linear_shape_function(xi, i);
2792 // see also BarycentricPolynomials<2>::compute_value
2793 case ReferenceCells::Triangle:
2794 {
2795 switch (i)
2796 {
2797 case 0:
2798 return 1.0 - xi[std::min(0, dim - 1)] -
2799 xi[std::min(1, dim - 1)];
2800 case 1:
2801 return xi[std::min(0, dim - 1)];
2802 case 2:
2803 return xi[std::min(1, dim - 1)];
2804 default:
2805 DEAL_II_ASSERT_UNREACHABLE();
2806 }
2807 }
2808 // see also BarycentricPolynomials<3>::compute_value
2809 case ReferenceCells::Tetrahedron:
2810 {
2811 switch (i)
2812 {
2813 case 0:
2814 return 1.0 - xi[std::min(0, dim - 1)] -
2815 xi[std::min(1, dim - 1)] - xi[std::min(2, dim - 1)];
2816 case 1:
2817 return xi[std::min(0, dim - 1)];
2818 case 2:
2819 return xi[std::min(1, dim - 1)];
2820 case 3:
2821 return xi[std::min(2, dim - 1)];
2822 default:
2823 DEAL_II_ASSERT_UNREACHABLE();
2824 }
2825 }
2826 // see also ScalarLagrangePolynomialPyramid::compute_value()
2827 case ReferenceCells::Pyramid:
2828 {
2829 const double Q14 = 0.25;
2830
2831 const double r = xi[std::min(0, dim - 1)];
2832 const double s = xi[std::min(1, dim - 1)];
2833 const double t = xi[std::min(2, dim - 1)];
2834
2835 const double ratio =
2836 (std::fabs(t - 1.0) > 1.0e-14 ? (r * s * t) / (1.0 - t) : 0.0);
2837
2838 if (i == 0)
2839 return Q14 * ((1.0 - r) * (1.0 - s) - t + ratio);
2840 if (i == 1)
2841 return Q14 * ((1.0 + r) * (1.0 - s) - t - ratio);
2842 if (i == 2)
2843 return Q14 * ((1.0 - r) * (1.0 + s) - t - ratio);
2844 if (i == 3)
2845 return Q14 * ((1.0 + r) * (1.0 + s) - t + ratio);
2846 else
2847 return t;
2848 }
2849 // see also ScalarLagrangePolynomialWedge::compute_value()
2850 case ReferenceCells::Wedge:
2851 return ReferenceCell(ReferenceCells::Triangle)
2852 .d_linear_shape_function<2>(Point<2>(xi[std::min(0, dim - 1)],
2853 xi[std::min(1, dim - 1)]),
2854 i % 3) *
2855 ReferenceCell(ReferenceCells::Line)
2856 .d_linear_shape_function<1>(Point<1>(xi[std::min(2, dim - 1)]),
2857 i / 3);
2858 default:
2859 DEAL_II_NOT_IMPLEMENTED();
2860 }
2861
2862 return 0.0;
2863}
2864
2865
2866
2867template <int dim>
2868inline Tensor<1, dim>
2869ReferenceCell::d_linear_shape_function_gradient(const Point<dim> &xi,
2870 const unsigned int i) const
2871{
2872 Assert(dim == get_dimension(),
2873 ExcMessage("You can only call this function with arguments for "
2874 "which 'dim' equals the dimension of the space in which "
2875 "the current reference cell lives. You have dim=" +
2876 std::to_string(dim) + " but the reference cell is " +
2877 std::to_string(get_dimension()) + " dimensional."));
2878 switch (this->kind)
2879 {
2880 case ReferenceCells::Vertex:
2881 case ReferenceCells::Line:
2882 case ReferenceCells::Quadrilateral:
2883 case ReferenceCells::Hexahedron:
2884 return GeometryInfo<dim>::d_linear_shape_function_gradient(xi, i);
2885 // see also BarycentricPolynomials<2>::compute_grad()
2886 case ReferenceCells::Triangle:
2887 switch (i)
2888 {
2889 case 0:
2890 return Point<dim>(-1.0, -1.0);
2891 case 1:
2892 return Point<dim>(+1.0, +0.0);
2893 case 2:
2894 return Point<dim>(+0.0, +1.0);
2895 default:
2896 DEAL_II_ASSERT_UNREACHABLE();
2897 }
2898 default:
2899 DEAL_II_NOT_IMPLEMENTED();
2900 }
2901
2902 return Point<dim>(+0.0, +0.0, +0.0);
2903}
2904
2905
2906
2907inline double
2908ReferenceCell::volume() const
2909{
2910 switch (this->kind)
2911 {
2912 case ReferenceCells::Vertex:
2913 return 0;
2914 case ReferenceCells::Line:
2915 return 1;
2916 case ReferenceCells::Triangle:
2917 return 1. / 2.;
2918 case ReferenceCells::Quadrilateral:
2919 return 1;
2920 case ReferenceCells::Tetrahedron:
2921 return 1. / 6.;
2922 case ReferenceCells::Pyramid:
2923 return 4. / 3.;
2924 case ReferenceCells::Wedge:
2925 return 1. / 2.;
2926 case ReferenceCells::Hexahedron:
2927 return 1;
2928 default:
2929 DEAL_II_NOT_IMPLEMENTED();
2930 }
2931
2932 return 0.0;
2933}
2934
2935
2936
2937template <int dim>
2938inline Point<dim>
2939ReferenceCell::barycenter() const
2940{
2941 Assert(dim == get_dimension(),
2942 ExcMessage("You can only call this function with arguments for "
2943 "which 'dim' equals the dimension of the space in which "
2944 "the current reference cell lives. You have dim=" +
2945 std::to_string(dim) + " but the reference cell is " +
2946 std::to_string(get_dimension()) + " dimensional."));
2947
2948 switch (this->kind)
2949 {
2950 case ReferenceCells::Vertex:
2951 return Point<dim>();
2952 case ReferenceCells::Line:
2953 return Point<dim>(1. / 2.);
2954 case ReferenceCells::Triangle:
2955 return Point<dim>(1. / 3., 1. / 3.);
2956 case ReferenceCells::Quadrilateral:
2957 return Point<dim>(1. / 2., 1. / 2.);
2958 case ReferenceCells::Tetrahedron:
2959 return Point<dim>(1. / 4., 1. / 4., 1. / 4.);
2960 case ReferenceCells::Pyramid:
2961 return Point<dim>(0, 0, 1. / 4.);
2962 case ReferenceCells::Wedge:
2963 return Point<dim>(1. / 3, 1. / 3, 1. / 2.);
2964 case ReferenceCells::Hexahedron:
2965 return Point<dim>(1. / 2., 1. / 2., 1. / 2.);
2966 default:
2967 DEAL_II_NOT_IMPLEMENTED();
2968 }
2969
2970 return Point<dim>();
2971}
2972
2973
2974
2975template <int dim>
2976inline bool
2977ReferenceCell::contains_point(const Point<dim> &p, const double tolerance) const
2978{
2979 Assert(dim == get_dimension(),
2980 ExcMessage("You can only call this function with arguments for "
2981 "which 'dim' equals the dimension of the space in which "
2982 "the current reference cell lives. You have dim=" +
2983 std::to_string(dim) + " but the reference cell is " +
2984 std::to_string(get_dimension()) + " dimensional."));
2985
2986 // Introduce abbreviations to silence compiler warnings about invalid
2987 // array accesses (that can't happen, of course, but the compiler
2988 // doesn't know that).
2989 constexpr unsigned int x_coordinate = 0;
2990 constexpr unsigned int y_coordinate = (dim >= 2 ? 1 : 0);
2991 constexpr unsigned int z_coordinate = (dim >= 3 ? 2 : 0);
2992
2993 switch (this->kind)
2994 {
2995 case ReferenceCells::Vertex:
2996 {
2997 // Vertices are special cases in that they do not actually
2998 // have coordinates. Error out if this function is called
2999 // with a vertex:
3000 Assert(false,
3001 ExcMessage("Vertices are zero-dimensional objects and "
3002 "as a consequence have no coordinates. You "
3003 "cannot meaningfully ask whether a point is "
3004 "inside a vertex (within a certain tolerance) "
3005 "without coordinate values."));
3006 return false;
3007 }
3008 case ReferenceCells::Line:
3009 case ReferenceCells::Quadrilateral:
3010 case ReferenceCells::Hexahedron:
3011 {
3012 for (unsigned int d = 0; d < dim; ++d)
3013 if ((p[d] < -tolerance) || (p[d] > 1 + tolerance))
3014 return false;
3015 return true;
3016 }
3017 case ReferenceCells::Triangle:
3018 case ReferenceCells::Tetrahedron:
3019 {
3020 // First make sure that we are in the first quadrant or octant
3021 for (unsigned int d = 0; d < dim; ++d)
3022 if (p[d] < -tolerance)
3023 return false;
3024
3025 // Now we also need to make sure that we are below the diagonal line
3026 // or plane that delineates the simplex. This diagonal is given by
3027 // sum(p[d])<=1, and a diagonal a distance eps away is given by
3028 // sum(p[d])<=1+eps*sqrt(d). (For example, the point at (1,1) is a
3029 // distance of 1/sqrt(2) away from the diagonal. That is, its
3030 // sum satisfies
3031 // sum(p[d]) = 2 <= 1 + (1/sqrt(2)) * sqrt(2)
3032 // in other words, it satisfies the predicate with eps=1/sqrt(2).)
3033 double sum = 0;
3034 for (unsigned int d = 0; d < dim; ++d)
3035 sum += p[d];
3036 return (sum <= 1 + tolerance * std::sqrt(1. * dim));
3037 }
3038 case ReferenceCells::Pyramid:
3039 {
3040 // A pyramid only lives in the upper half-space:
3041 if (p[z_coordinate] < -tolerance)
3042 return false;
3043
3044 // It also only lives in the space below z=1:
3045 if (p[z_coordinate] > 1 + tolerance)
3046 return false;
3047
3048 // Within what's left of the space, a pyramid is a cone that tapers
3049 // towards the top. First compute the distance of the point to the
3050 // axis in the max norm (this is the right norm because the vertices
3051 // of the pyramid are at points +/-1, +/-1):
3052 const double distance_from_axis =
3053 std::max(std::fabs(p[x_coordinate]), std::fabs(p[y_coordinate]));
3054
3055 // We are inside the pyramid if the distance from the axis is less
3056 // than (1-z)
3057 return (distance_from_axis <= 1 + tolerance - p[z_coordinate]);
3058 }
3059 case ReferenceCells::Wedge:
3060 {
3061 // The wedge we use is a triangle extruded into the third
3062 // dimension by one unit. So we can use the same logic as for
3063 // triangles above (i.e., for the simplex above, using dim==2)
3064 // and then check the third dimension separately.
3065
3066 if ((p[x_coordinate] < -tolerance) || (p[y_coordinate] < -tolerance))
3067 return false;
3068
3069 const double sum = p[x_coordinate] + p[y_coordinate];
3070 if (sum > 1 + tolerance * std::sqrt(2.0))
3071 return false;
3072
3073 if (p[z_coordinate] < -tolerance)
3074 return false;
3075 if (p[z_coordinate] > 1 + tolerance)
3076 return false;
3077
3078 return true;
3079 }
3080 default:
3081 DEAL_II_NOT_IMPLEMENTED();
3082 }
3083
3084 return false;
3085}
3086
3087
3088
3089template <int dim>
3090inline Tensor<1, dim>
3091ReferenceCell::unit_tangential_vectors(const unsigned int face_no,
3092 const unsigned int i) const
3093{
3094 Assert(dim == get_dimension(),
3095 ExcMessage("You can only call this function with arguments for "
3096 "which 'dim' equals the dimension of the space in which "
3097 "the current reference cell lives. You have dim=" +
3098 std::to_string(dim) + " but the reference cell is " +
3099 std::to_string(get_dimension()) + " dimensional."));
3100 AssertIndexRange(i, dim - 1);
3101
3102 switch (this->kind)
3103 {
3104 case ReferenceCells::Vertex:
3105 case ReferenceCells::Line:
3106 case ReferenceCells::Quadrilateral:
3107 case ReferenceCells::Hexahedron:
3108 AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
3109 return GeometryInfo<dim>::unit_tangential_vectors[face_no][i];
3110 case ReferenceCells::Triangle:
3111 {
3112 AssertIndexRange(face_no, 3);
3113 static const std::array<Tensor<1, dim>, 3> table = {
3114 {Point<dim>(1, 0),
3115 Point<dim>(-std::sqrt(0.5), +std::sqrt(0.5)),
3116 Point<dim>(0, -1)}};
3117
3118 return table[face_no];
3119 }
3120 case ReferenceCells::Tetrahedron:
3121 {
3122 AssertIndexRange(face_no, 4);
3123 static const ndarray<Tensor<1, dim>, 4, 2> table = {
3124 {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
3125 {{Point<dim>(1, 0, 0), Point<dim>(0, 0, 1)}},
3126 {{Point<dim>(0, 0, 1), Point<dim>(0, 1, 0)}},
3127 {{Point<dim>(-std::pow(1.0 / 3.0, 1.0 / 4.0),
3128 +std::pow(1.0 / 3.0, 1.0 / 4.0),
3129 0),
3130 Point<dim>(-std::pow(1.0 / 3.0, 1.0 / 4.0),
3131 0,
3132 +std::pow(1.0 / 3.0, 1.0 / 4.0))}}}};
3133
3134 return table[face_no][i];
3135 }
3136 case ReferenceCells::Pyramid:
3137 {
3138 AssertIndexRange(face_no, 5);
3139 static const ndarray<Tensor<1, dim>, 5, 2> table = {
3140 {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
3141 {{Point<dim>(+1.0 / sqrt(2.0), 0, +1.0 / sqrt(2.0)),
3142 Point<dim>(0, 1, 0)}},
3143 {{Point<dim>(+1.0 / sqrt(2.0), 0, -1.0 / sqrt(2.0)),
3144 Point<dim>(0, 1, 0)}},
3145 {{Point<dim>(1, 0, 0),
3146 Point<dim>(0, +1.0 / sqrt(2.0), +1.0 / sqrt(2.0))}},
3147 {{Point<dim>(1, 0, 0),
3148 Point<dim>(0, +1.0 / sqrt(2.0), -1.0 / sqrt(2.0))}}}};
3149
3150 return table[face_no][i];
3151 }
3152 case ReferenceCells::Wedge:
3153 {
3154 AssertIndexRange(face_no, 5);
3155 static const ndarray<Tensor<1, dim>, 5, 2> table = {
3156 {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
3157 {{Point<dim>(1, 0, 0), Point<dim>(0, 1, 0)}},
3158 {{Point<dim>(1, 0, 0), Point<dim>(0, 0, 1)}},
3159 {{Point<dim>(-1 / std::sqrt(2.0), +1 / std::sqrt(2.0), 0),
3160 Point<dim>(0, 0, 1)}},
3161 {{Point<dim>(0, 0, 1), Point<dim>(0, 1, 0)}}}};
3162
3163 return table[face_no][i];
3164 }
3165 default:
3166 DEAL_II_NOT_IMPLEMENTED();
3167 }
3168
3169
3170 return {};
3171}
3172
3173
3174
3175template <int dim>
3176inline Tensor<1, dim>
3177ReferenceCell::unit_normal_vectors(const unsigned int face_no) const
3178{
3179 Assert(dim == get_dimension(),
3180 ExcMessage("You can only call this function with arguments for "
3181 "which 'dim' equals the dimension of the space in which "
3182 "the current reference cell lives. You have dim=" +
3183 std::to_string(dim) + " but the reference cell is " +
3184 std::to_string(get_dimension()) + " dimensional."));
3185
3186 if (is_hyper_cube())
3187 {
3188 AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
3189 return GeometryInfo<dim>::unit_normal_vector[face_no];
3190 }
3191 else if (dim == 2)
3192 {
3193 Assert(*this == ReferenceCells::Triangle, ExcInternalError());
3194
3195 // Return the rotated vector
3196 return cross_product_2d(unit_tangential_vectors<dim>(face_no, 0));
3197 }
3198 else if (dim == 3)
3199 {
3200 return cross_product_3d(unit_tangential_vectors<dim>(face_no, 0),
3201 unit_tangential_vectors<dim>(face_no, 1));
3202 }
3203
3204 DEAL_II_NOT_IMPLEMENTED();
3205
3206 return {};
3207}
3208
3209
3210
3211inline unsigned int
3212ReferenceCell::n_face_orientations(const unsigned int face_no) const
3213{
3214 AssertIndexRange(face_no, n_faces());
3215 if (get_dimension() == 1)
3216 return 1;
3217 if (get_dimension() == 2)
3218 return 2;
3219 else if (face_reference_cell(face_no) == ReferenceCells::Quadrilateral)
3220 return 8;
3221 else if (face_reference_cell(face_no) == ReferenceCells::Triangle)
3222 return 6;
3223
3224 DEAL_II_ASSERT_UNREACHABLE();
3225 return numbers::invalid_unsigned_int;
3226}
3227
3228
3229
3230inline bool
3231ReferenceCell::standard_vs_true_line_orientation(
3232 const unsigned int line,
3233 const unsigned int face,
3234 const types::geometric_orientation combined_face_orientation,
3235 const bool line_orientation) const
3236{
3237 if (*this == ReferenceCells::Hexahedron)
3238 {
3239 static constexpr ::ndarray<bool, 2, 8> bool_table{
3240 {{{true, true, false, true, false, false, true, false}},
3241 {{true, true, true, false, false, false, false, true}}}};
3242
3243 return (line_orientation ==
3244 bool_table[line / 2][combined_face_orientation]);
3245 }
3246 else if (*this == ReferenceCells::Tetrahedron)
3247 {
3248 static constexpr unsigned int X = numbers::invalid_unsigned_int;
3249 static constexpr ::ndarray<unsigned int, 4, 3> combined_lines{
3250 {{{0, 0, 0}}, {{X, 0, 1}}, {{X, 0, X}}, {{X, X, X}}}};
3251
3252 const auto combined_line = combined_lines[face][line];
3253
3254 Assert(combined_line != X,
3255 ExcMessage(
3256 "This function can only be called for following face-line "
3257 "combinations: (0,0), (0,1), (0,2), (1,1), (1,2), (2,1),"));
3258
3259 static constexpr ::ndarray<bool, 2, 6> bool_table{
3260 {{{false, true, false, true, false, true}},
3261 {{true, false, true, false, true, false}}}};
3262
3263 return (line_orientation ==
3264 bool_table[combined_line][combined_face_orientation]);
3265 }
3266 else
3267 // TODO: This might actually be wrong for some of the other
3268 // kinds of objects. We should check this
3269 return true;
3270}
3271
3272
3273
3274namespace internal
3275{
3276 template <typename T>
3277 class NoPermutation : public ::ExceptionBase
3278 {
3279 public:
3293
3297 virtual ~NoPermutation() noexcept override = default;
3298
3302 virtual void
3303 print_info(std::ostream &out) const override
3304 {
3305 out << '[';
3306
3307 const unsigned int n_vertices = entity_type.n_vertices();
3308
3309 for (unsigned int i = 0; i < n_vertices; ++i)
3310 {
3311 out << vertices_0[i];
3312 if (i + 1 != n_vertices)
3313 out << ',';
3314 }
3315
3316 out << "] is not a valid permutation of [";
3317
3318 for (unsigned int i = 0; i < n_vertices; ++i)
3319 {
3320 out << vertices_1[i];
3321 if (i + 1 != n_vertices)
3322 out << ',';
3323 }
3324
3325 out << "]." << std::endl;
3326 }
3327
3331 const ReferenceCell entity_type;
3332
3336 const ArrayView<const T> vertices_0;
3337
3341 const ArrayView<const T> vertices_1;
3342 };
3343
3351 DeclException2(
3352 ExcNonMatchingReferenceCellTypes,
3353 ReferenceCell,
3354 ReferenceCell,
3355 << "The reference-cell type used on this cell (" << arg1.to_string()
3356 << ") does not match the reference-cell type of the finite element "
3357 << "associated with this cell (" << arg2.to_string() << "). "
3358 << "Did you accidentally use simplex elements on hypercube meshes "
3359 << "(or the other way around), or are you using a mixed mesh and "
3360 << "assigned a simplex element to a hypercube cell (or the other "
3361 << "way around) via the active_fe_index?");
3362} // namespace internal
3363
3364
3365
3366template <typename T, std::size_t N>
3367inline types::geometric_orientation
3368ReferenceCell::compute_orientation(const std::array<T, N> &vertices_0,
3369 const std::array<T, N> &vertices_1) const
3370{
3371 Assert(N >= n_vertices(),
3372 ExcMessage("The number of array elements must be equal to or "
3373 "greater than the number of vertices of the cell "
3374 "referenced by this object."));
3375
3376 // Call the non-deprecated function, taking care of calling it only with
3377 // those array elements that we actually care about (see the note
3378 // in the documentation about the arguments potentially being
3379 // larger arrays than necessary).
3380 return get_combined_orientation(
3381 make_array_view(vertices_0.begin(), vertices_0.begin() + n_vertices()),
3382 make_array_view(vertices_1.begin(), vertices_1.begin() + n_vertices()));
3383}
3384
3385
3386
3387template <typename T>
3388types::geometric_orientation
3389ReferenceCell::get_combined_orientation(
3390 const ArrayView<const T> &vertices_0,
3391 const ArrayView<const T> &vertices_1) const
3392{
3393 Assert(vertices_0.size() == n_vertices(),
3394 ExcMessage("The number of array elements must be equal to "
3395 "the number of vertices of the cell "
3396 "referenced by this object."));
3397 Assert(vertices_1.size() == n_vertices(),
3398 ExcMessage("The number of array elements must be equal to "
3399 "the number of vertices of the cell "
3400 "referenced by this object."));
3401
3402 auto compute_orientation = [&](const auto &table) {
3403 for (types::geometric_orientation o = 0; o < table.size(); ++o)
3404 {
3405 bool match = true;
3406 for (unsigned int j = 0; j < table[o].size(); ++j)
3407 match = (match && vertices_0[j] == vertices_1[table[o][j]]);
3408
3409 if (match)
3410 return o;
3411 }
3412
3413 Assert(false, (internal::NoPermutation<T>(*this, vertices_0, vertices_1)));
3414 return std::numeric_limits<types::geometric_orientation>::max();
3415 };
3416
3417 switch (this->kind)
3418 {
3419 case ReferenceCells::Vertex:
3420 // TODO: we can get rid of this special-case and use
3421 // vertex_vertex_permutations once we make 0 the default orientation.
3422 (void)vertex_vertex_permutations;
3423 if (vertices_0[0] == vertices_1[0])
3424 return default_combined_face_orientation();
3425 break;
3426 case ReferenceCells::Line:
3427 return compute_orientation(line_vertex_permutations);
3428 case ReferenceCells::Triangle:
3429 return compute_orientation(triangle_vertex_permutations);
3430 case ReferenceCells::Quadrilateral:
3431 return compute_orientation(quadrilateral_vertex_permutations);
3432 default:
3433 DEAL_II_NOT_IMPLEMENTED();
3434 }
3435
3436 Assert(false, (internal::NoPermutation<T>(*this, vertices_0, vertices_1)));
3437 return std::numeric_limits<types::geometric_orientation>::max();
3438}
3439
3440
3441
3442template <typename T, std::size_t N>
3443inline std::array<T, N>
3444ReferenceCell::permute_according_orientation(
3445 const std::array<T, N> &vertices,
3446 const unsigned int orientation) const
3447{
3448 Assert(N >= n_vertices(),
3449 ExcMessage("The number of array elements must be equal to or "
3450 "greater than the number of vertices of the cell "
3451 "referenced by this object."));
3452
3453 // Call the non-deprecated function, taking care of calling it only with
3454 // those array elements that we actually care about (see the note
3455 // in the documentation about the arguments potentially being
3456 // larger arrays than necessary).
3457 const auto permutation = permute_by_combined_orientation(
3458 make_array_view(vertices.begin(), vertices.begin() + n_vertices()),
3459 orientation);
3460
3461 std::array<T, N> temp;
3462 std::copy(permutation.begin(), permutation.end(), temp.begin());
3463
3464 return temp;
3465}
3466
3467
3468
3469template <typename T>
3470boost::container::small_vector<T, 8>
3471ReferenceCell::permute_by_combined_orientation(
3472 const ArrayView<const T> &vertices,
3473 const types::geometric_orientation orientation) const
3474{
3475 AssertDimension(vertices.size(), n_vertices());
3476 boost::container::small_vector<T, 8> result(n_vertices());
3477
3478 auto permute = [&](const auto &table) {
3479 AssertIndexRange(orientation, table.size());
3480 for (unsigned int j = 0; j < table[orientation].size(); ++j)
3481 result[j] = vertices[table[orientation][j]];
3482 };
3483
3484 switch (this->kind)
3485 {
3486 case ReferenceCells::Vertex:
3487 // TODO: we can get rid of this special-case and use
3488 // vertex_vertex_permutations once we make 0 the default orientation.
3489 std::copy(vertices.begin(), vertices.end(), result.begin());
3490 break;
3491 case ReferenceCells::Line:
3492 permute(line_vertex_permutations);
3493 break;
3494 case ReferenceCells::Triangle:
3495 permute(triangle_vertex_permutations);
3496 break;
3497 case ReferenceCells::Quadrilateral:
3498 permute(quadrilateral_vertex_permutations);
3499 break;
3500 default:
3501 DEAL_II_NOT_IMPLEMENTED();
3502 }
3503
3504 return result;
3505}
3506
3507
3508
3509inline types::geometric_orientation
3510ReferenceCell::get_inverse_combined_orientation(
3511 const types::geometric_orientation orientation) const
3512{
3513 switch (this->kind)
3514 {
3515 case ReferenceCells::Vertex:
3516 // Things are always default-oriented in 1D
3517 return orientation;
3518
3519 case ReferenceCells::Line:
3520 // the 1d orientations are the identity and a flip: i.e., the identity
3521 // and an involutory mapping
3522 return orientation;
3523
3524 case ReferenceCells::Triangle:
3525 {
3526 AssertIndexRange(orientation, 6);
3527 constexpr std::array<types::geometric_orientation, 6> inverses{
3528 {0, 1, 2, 5, 4, 3}};
3529 return inverses[orientation];
3530 }
3531 case ReferenceCells::Quadrilateral:
3532 {
3533 AssertIndexRange(orientation, 8);
3534 constexpr std::array<types::geometric_orientation, 8> inverses{
3535 {0, 1, 2, 7, 4, 5, 6, 3}};
3536 return inverses[orientation];
3537 }
3538 default:
3539 DEAL_II_NOT_IMPLEMENTED();
3540 }
3541
3542 return std::numeric_limits<types::geometric_orientation>::max();
3543}
3544
3545
3546
3547template <>
3548unsigned int
3549ReferenceCell::vtk_lexicographic_to_node_index<0>(
3550 const std::array<unsigned, 0> &node_indices,
3551 const std::array<unsigned, 0> &nodes_per_direction,
3552 const bool legacy_format) const;
3553
3554template <>
3555unsigned int
3556ReferenceCell::vtk_lexicographic_to_node_index<1>(
3557 const std::array<unsigned, 1> &node_indices,
3558 const std::array<unsigned, 1> &nodes_per_direction,
3559 const bool legacy_format) const;
3560
3561template <>
3562unsigned int
3563ReferenceCell::vtk_lexicographic_to_node_index<2>(
3564 const std::array<unsigned, 2> &node_indices,
3565 const std::array<unsigned, 2> &nodes_per_direction,
3566 const bool legacy_format) const;
3567
3568template <>
3569unsigned int
3570ReferenceCell::vtk_lexicographic_to_node_index<3>(
3571 const std::array<unsigned, 3> &node_indices,
3572 const std::array<unsigned, 3> &nodes_per_direction,
3573 const bool legacy_format) const;
3574
3575DEAL_II_NAMESPACE_CLOSE
3576
3577#endif
array_view.h
make_array_view
ArrayView< std::remove_reference_t< typename std::iterator_traits< Iterator >::reference >, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Definition array_view.h:949
ArrayView::begin
iterator begin() const
Definition array_view.h:702
ArrayView::end
iterator end() const
Definition array_view.h:711
ArrayView::size
std::size_t size() const
Definition array_view.h:684
Mapping
Abstract base class for mapping classes.
Definition mapping.h:318
Point
Definition point.h:111
Point::unit_vector
static constexpr Point< dim, Number > unit_vector(const unsigned int i)
ReferenceCell::faces_for_given_vertex
ArrayView< const unsigned int > faces_for_given_vertex(const unsigned int vertex_index) const
Definition reference_cell.h:1250
ReferenceCell::vertex_indices
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
Definition reference_cell.h:1984
ReferenceCell::reversed_combined_line_orientation
static constexpr types::geometric_orientation reversed_combined_line_orientation()
Definition reference_cell.h:2048
ReferenceCell::n_vertices
unsigned int n_vertices() const
Definition reference_cell.h:1645
ReferenceCell::get_gauss_type_quadrature
Quadrature< dim > get_gauss_type_quadrature(const unsigned n_points_1d) const
Definition reference_cell.cc:186
ReferenceCell::is_hyper_cube
bool is_hyper_cube() const
Definition reference_cell.h:1582
ReferenceCell::vtk_linear_type
unsigned int vtk_linear_type() const
Definition reference_cell.cc:379
ReferenceCell::serialize
void serialize(Archive &archive, const unsigned int)
Definition reference_cell.h:1234
ReferenceCell::n_isotropic_refinement_choices
unsigned int n_isotropic_refinement_choices() const
Definition reference_cell.h:1868
ReferenceCell::standard_to_real_face_vertex
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const types::geometric_orientation face_orientation) const
Definition reference_cell.h:2535
ReferenceCell::compute_orientation
types::geometric_orientation compute_orientation(const std::array< T, N > &vertices_0, const std::array< T, N > &vertices_1) const
Definition reference_cell.h:3368
ReferenceCell::standard_vertex_to_face_and_vertex_index
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const
Definition reference_cell.h:2124
ReferenceCell::vertex
Point< dim > vertex(const unsigned int v) const
Definition reference_cell.h:1706
ReferenceCell::get_nodal_type_quadrature
const Quadrature< dim > & get_nodal_type_quadrature() const
Definition reference_cell.cc:208
ReferenceCell::face_vertex_location
Point< dim > face_vertex_location(const unsigned int face, const unsigned int vertex) const
Definition reference_cell.h:2519
ReferenceCell::new_isotropic_child_face_lines
constexpr ::ndarray< unsigned int, 12, 4 > new_isotropic_child_face_lines(const unsigned int refinement_choice) const
Definition reference_cell.h:1307
ReferenceCell::face_to_cell_lines
unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const types::geometric_orientation face_orientation) const
Definition reference_cell.h:2330
ReferenceCell::vertex_vertex_permutations
static constexpr ndarray< unsigned int, 1, 1 > vertex_vertex_permutations
Definition reference_cell.h:1017
ReferenceCell::isotropic_child_indices
std_cxx20::ranges::iota_view< unsigned int, unsigned int > isotropic_child_indices() const
Definition reference_cell.h:1976
ReferenceCell::standard_line_to_face_and_line_index
std::array< unsigned int, 2 > standard_line_to_face_and_line_index(const unsigned int line) const
Definition reference_cell.h:2186
ReferenceCell::new_isotropic_child_face_line_vertices
constexpr ::ndarray< unsigned int, 12, 4, 2 > new_isotropic_child_face_line_vertices(const unsigned int refinement_choice) const
Definition reference_cell.h:1386
ReferenceCell::d_linear_shape_function
double d_linear_shape_function(const Point< dim > &xi, const unsigned int i) const
Definition reference_cell.h:2775
ReferenceCell::get_inverse_combined_orientation
types::geometric_orientation get_inverse_combined_orientation(const types::geometric_orientation orientation) const
Definition reference_cell.h:3510
ReferenceCell::unit_tangential_vectors
Tensor< 1, dim > unit_tangential_vectors(const unsigned int face_no, const unsigned int i) const
Definition reference_cell.h:3091
ReferenceCell::line_vertex_permutations
static constexpr ndarray< unsigned int, 2, 2 > line_vertex_permutations
Definition reference_cell.h:1023
ReferenceCell::vtk_lexicographic_to_node_index
unsigned int vtk_lexicographic_to_node_index(const std::array< unsigned, dim > &node_indices, const std::array< unsigned, dim > &nodes_per_direction, const bool legacy_format) const
ReferenceCell::permute_according_orientation
std::array< T, N > permute_according_orientation(const std::array< T, N > &vertices, const unsigned int orientation) const
Definition reference_cell.h:3444
ReferenceCell::new_isotropic_child_cell_vertices
constexpr ::ndarray< unsigned int, 8, 4 > new_isotropic_child_cell_vertices(const unsigned int refinement_choice) const
Definition reference_cell.h:1533
ReferenceCell::contains_point
bool contains_point(const Point< dim > &p, const double tolerance=0) const
Definition reference_cell.h:2977
ReferenceCell::new_isotropic_child_cell_faces
constexpr ::ndarray< unsigned int, 8, 6 > new_isotropic_child_cell_faces(const unsigned int refinement_choice) const
Definition reference_cell.h:1465
ReferenceCell::n_isotropic_children
unsigned int n_isotropic_children() const
Definition reference_cell.h:1924
ReferenceCell::permute_by_combined_orientation
boost::container::small_vector< T, 8 > permute_by_combined_orientation(const ArrayView< const T > &vertices, const types::geometric_orientation orientation) const
Definition reference_cell.h:3471
ReferenceCell::gmsh_element_type
unsigned int gmsh_element_type() const
Definition reference_cell.cc:727
ReferenceCell::n_face_orientations
unsigned int n_face_orientations(const unsigned int face_no) const
Definition reference_cell.h:3212
ReferenceCell::quadrilateral_vertex_permutations
static constexpr ndarray< unsigned int, 8, 4 > quadrilateral_vertex_permutations
Definition reference_cell.h:1041
ReferenceCell::volume
double volume() const
Definition reference_cell.h:2908
ReferenceCell::n_faces
unsigned int n_faces() const
Definition reference_cell.h:1830
ReferenceCell::line_to_cell_vertices
unsigned int line_to_cell_vertices(const unsigned int line, const unsigned int vertex) const
Definition reference_cell.h:2248
ReferenceCell::kind
std::uint8_t kind
Definition reference_cell.h:1002
ReferenceCell::ReferenceCell
constexpr ReferenceCell()
Definition reference_cell.h:1226
ReferenceCell::exodusii_vertex_to_deal_vertex
unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const
Definition reference_cell.cc:253
ReferenceCell::standard_vs_true_line_orientation
bool standard_vs_true_line_orientation(const unsigned int line, const unsigned int face, const types::geometric_orientation face_orientation, const bool line_orientation) const
Definition reference_cell.h:3231
ReferenceCell::unv_vertex_to_deal_vertex
unsigned int unv_vertex_to_deal_vertex(const unsigned int vertex_n) const
Definition reference_cell.cc:343
ReferenceCell::default_combined_face_orientation
static constexpr types::geometric_orientation default_combined_face_orientation()
Definition reference_cell.h:2037
ReferenceCell::vtk_vertex_to_deal_vertex
unsigned int vtk_vertex_to_deal_vertex(const unsigned int vertex_index) const
Definition reference_cell.cc:672
ReferenceCell::get_combined_orientation
types::geometric_orientation get_combined_orientation(const ArrayView< const T > &vertices_0, const ArrayView< const T > &vertices_1) const
Definition reference_cell.h:3389
ReferenceCell::get_isotropic_refinement_choice
IsotropicRefinementChoice get_isotropic_refinement_choice(const unsigned int ref_choice) const
Definition reference_cell.h:1898
ReferenceCell::face_to_cell_vertices
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const types::geometric_orientation face_orientation) const
Definition reference_cell.h:2417
ReferenceCell::get_default_mapping
std::unique_ptr< Mapping< dim, spacedim > > get_default_mapping(const unsigned int degree) const
Definition reference_cell.cc:121
ReferenceCell::get_midpoint_quadrature
Quadrature< dim > get_midpoint_quadrature() const
Definition reference_cell.h:1629
ReferenceCell::vtk_quadratic_type
unsigned int vtk_quadratic_type() const
Definition reference_cell.cc:411
ReferenceCell::n_lines
unsigned int n_lines() const
Definition reference_cell.h:1675
ReferenceCell::vtk_lagrange_type
unsigned int vtk_lagrange_type() const
Definition reference_cell.cc:443
ReferenceCell::line_indices
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
Definition reference_cell.h:1992
ReferenceCell::get_dimension
unsigned int get_dimension() const
Definition reference_cell.h:1602
ReferenceCell::n_children
unsigned int n_children(const RefinementCase< dim > ref_case=RefinementCase< dim >::isotropic_refinement) const
Definition reference_cell.h:1957
ReferenceCell::face_reference_cell
ReferenceCell face_reference_cell(const unsigned int face_no) const
Definition reference_cell.h:2000
ReferenceCell::exodusii_face_to_deal_face
unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const
Definition reference_cell.cc:296
ReferenceCell::n_vertices_to_type
static ReferenceCell n_vertices_to_type(const int dim, const unsigned int n_vertices)
Definition reference_cell.h:2733
ReferenceCell::memory_consumption
static constexpr std::size_t memory_consumption()
Definition reference_cell.h:1242
ReferenceCell::operator>>
friend std::istream & operator>>(std::istream &in, ReferenceCell &reference_cell)
Definition reference_cell.cc:1145
ReferenceCell::get_default_linear_mapping
const Mapping< dim, spacedim > & get_default_linear_mapping() const
Definition reference_cell.cc:148
ReferenceCell::to_string
std::string to_string() const
Definition reference_cell.cc:88
ReferenceCell::is_simplex
bool is_simplex() const
Definition reference_cell.h:1592
ReferenceCell::child_cell_on_face
unsigned int child_cell_on_face(const unsigned int face, const unsigned int subface, const types::geometric_orientation face_orientation=default_combined_face_orientation()) const
Definition reference_cell.h:2058
ReferenceCell::operator!=
constexpr bool operator!=(const ReferenceCell &type) const
Definition reference_cell.h:1108
ReferenceCell::operator==
constexpr bool operator==(const ReferenceCell &type) const
Definition reference_cell.h:1100
ReferenceCell::barycenter
Point< dim > barycenter() const
Definition reference_cell.h:2939
ReferenceCell::standard_to_real_face_line
unsigned int standard_to_real_face_line(const unsigned int line, const unsigned int face, const types::geometric_orientation face_orientation) const
Definition reference_cell.h:2583
ReferenceCell::unit_normal_vectors
Tensor< 1, dim > unit_normal_vectors(const unsigned int face_no) const
Definition reference_cell.h:3177
ReferenceCell::face_indices
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
Definition reference_cell.h:1860
ReferenceCell::operator<<
friend std::ostream & operator<<(std::ostream &out, const ReferenceCell &reference_cell)
Definition reference_cell.cc:1132
ReferenceCell::closest_point
Point< dim > closest_point(const Point< dim > &p) const
Definition reference_cell.cc:950
ReferenceCell::triangle_vertex_permutations
static constexpr ndarray< unsigned int, 6, 3 > triangle_vertex_permutations
Definition reference_cell.h:1029
ReferenceCell::d_linear_shape_function_gradient
Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i) const
Definition reference_cell.h:2869
internal::NoPermutation::print_info
virtual void print_info(std::ostream &out) const override
Definition reference_cell.h:3303
internal::NoPermutation::entity_type
const ReferenceCell entity_type
Definition reference_cell.h:3331
internal::NoPermutation::~NoPermutation
virtual ~NoPermutation() noexcept override=default
internal::NoPermutation::vertices_0
const ArrayView< const T > vertices_0
Definition reference_cell.h:3336
internal::NoPermutation::vertices_1
const ArrayView< const T > vertices_1
Definition reference_cell.h:3341
internal::NoPermutation::NoPermutation
NoPermutation(const ReferenceCell &entity_type, const ArrayView< const T > &vertices_0, const ArrayView< const T > &vertices_1)
Definition reference_cell.h:3283
DEAL_II_DEPRECATED
#define DEAL_II_DEPRECATED
Definition config.h:205
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
geometry_info.h
internal::ExcNonMatchingReferenceCellTypes
static ::ExceptionBase & ExcNonMatchingReferenceCellTypes(ReferenceCell arg1, ReferenceCell arg2)
Assert
#define Assert(cond, exc)
Definition exceptions.h:1659
DeclException2
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:534
AssertDimension
#define AssertDimension(dim1, dim2)
Definition exceptions.h:2026
AssertIndexRange
#define AssertIndexRange(index, range)
Definition exceptions.h:2115
StandardExceptions::ExcInternalError
static ::ExceptionBase & ExcInternalError()
StandardExceptions::ExcMessage
static ::ExceptionBase & ExcMessage(std::string arg1)
DEAL_II_ASSERT_UNREACHABLE
#define DEAL_II_ASSERT_UNREACHABLE()
Definition exceptions.h:1938
DEAL_II_NOT_IMPLEMENTED
#define DEAL_II_NOT_IMPLEMENTED()
Definition exceptions.h:1855
ReferenceCells::Tetrahedron
constexpr const ReferenceCell Tetrahedron
Definition reference_cell.h:1152
ReferenceCells::Quadrilateral
constexpr const ReferenceCell Quadrilateral
Definition reference_cell.h:1150
ReferenceCells::Wedge
constexpr const ReferenceCell Wedge
Definition reference_cell.h:1156
ReferenceCells::Pyramid
constexpr const ReferenceCell Pyramid
Definition reference_cell.h:1154
ReferenceCells::Invalid
constexpr const ReferenceCell Invalid
Definition reference_cell.h:1160
ReferenceCells::Triangle
constexpr const ReferenceCell Triangle
Definition reference_cell.h:1148
ReferenceCells::Hexahedron
constexpr const ReferenceCell Hexahedron
Definition reference_cell.h:1158
ReferenceCells::max_n_lines
constexpr unsigned int max_n_lines()
Definition reference_cell.h:2716
ReferenceCells::get_hypercube
constexpr const ReferenceCell & get_hypercube()
Definition reference_cell.h:2685
ReferenceCells::Vertex
constexpr const ReferenceCell Vertex
Definition reference_cell.h:1144
ReferenceCells::max_n_vertices
constexpr unsigned int max_n_vertices()
Definition reference_cell.h:2707
ReferenceCells::max_n_faces
constexpr unsigned int max_n_faces()
Definition reference_cell.h:2725
ReferenceCells::Line
constexpr const ReferenceCell Line
Definition reference_cell.h:1146
ReferenceCells::get_simplex
constexpr const ReferenceCell & get_simplex()
Definition reference_cell.h:2663
internal::make_reference_cell_from_int
constexpr ReferenceCell make_reference_cell_from_int(const std::uint8_t kind)
Definition reference_cell.h:1118
internal::split_face_orientation
std::tuple< bool, bool, bool > split_face_orientation(const types::geometric_orientation combined_face_orientation)
Definition tria_orientation.h:45
numbers::invalid_unsigned_int
static const unsigned int invalid_unsigned_int
Definition types.h:232
std_cxx20::ranges::iota_view
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:45
std
STL namespace.
std::min
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
Definition vectorization.h:6960
std::max
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
Definition vectorization.h:6943
std::sqrt
::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &)
Definition vectorization.h:6869
std::pow
::VectorizedArray< Number, width > pow(const ::VectorizedArray< Number, width > &, const Number p)
Definition vectorization.h:6885
types::geometric_orientation
unsigned char geometric_orientation
Definition types.h:40
ndarray
typename internal::ndarray::HelperArray< T, Ns... >::type ndarray
Definition ndarray.h:107
IsotropicRefinementChoice
IsotropicRefinementChoice
Definition reference_cell.h:75
IsotropicRefinementChoice::cut_tet_49
@ cut_tet_49
IsotropicRefinementChoice::cut_tet_57
@ cut_tet_57
IsotropicRefinementChoice::cut_tet_68
@ cut_tet_68
IsotropicRefinementChoice::isotropic_refinement
@ isotropic_refinement
operator>>
std::istream & operator>>(std::istream &in, ReferenceCell &reference_cell)
Definition reference_cell.cc:1145
operator<<
std::ostream & operator<<(std::ostream &out, const ReferenceCell &reference_cell)
Definition reference_cell.cc:1132
GeometryInfo::child_cell_on_face
static unsigned int child_cell_on_face(const RefinementCase< dim > &ref_case, const unsigned int face, const unsigned int subface, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement)
GeometryInfo::face_to_cell_vertices
static unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
GeometryInfo::standard_hex_line_to_quad_line_index
static std::array< unsigned int, 2 > standard_hex_line_to_quad_line_index(const unsigned int line)
GeometryInfo::standard_to_real_face_line
static unsigned int standard_to_real_face_line(const unsigned int line, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
GeometryInfo::standard_hex_vertex_to_quad_vertex_index
static std::array< unsigned int, 2 > standard_hex_vertex_to_quad_vertex_index(const unsigned int vertex)
GeometryInfo::standard_quad_vertex_to_line_vertex_index
static std::array< unsigned int, 2 > standard_quad_vertex_to_line_vertex_index(const unsigned int vertex)
GeometryInfo::d_linear_shape_function
static double d_linear_shape_function(const Point< dim > &xi, const unsigned int i)
GeometryInfo::n_children
static unsigned int n_children(const RefinementCase< dim > &refinement_case)
GeometryInfo::face_to_cell_lines
static unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
GeometryInfo::d_linear_shape_function_gradient
static Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i)
tria_orientation.h
Generated by doxygen 1.9.8