deal.II version GIT relicensing-6809-ge913b9bb34 2026-09-25 17:20:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
reference_cell.h
Go to the documentation of this file.
1// -----------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception OR LGPL-2.1-or-later
4// Copyright (C) 2020 - 2026 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Detailed license information governing the source code and contributions
9// can be found in LICENSE.md and CONTRIBUTING.md at the top level directory.
10//
11// -----------------------------------------------------------------------------
12
13#ifndef dealii_tria_reference_cell_h
14#define dealii_tria_reference_cell_h
15
16#include <deal.II/base/config.h>
17
18#include <deal.II/base/array_view.h>
19#include <deal.II/base/geometry_info.h>
20#include <deal.II/base/ndarray.h>
21#include <deal.II/base/point.h>
22#include <deal.II/base/tensor.h>
23#include <deal.II/base/types.h>
24#include <deal.II/base/utilities.h>
25
26#include <deal.II/grid/tria_orientation.h>
27
28#include <boost/container/small_vector.hpp>
29
30#include <iosfwd>
31#include <string>
32
33DEAL_II_NAMESPACE_OPEN
34
35// Forward declarations
36#ifndef DOXYGEN
37template <int dim, int spacedim>
38class Mapping;
39
40template <int dim>
41class Quadrature;
42
43template <int dim>
44class ReferenceCell;
45#endif
46
47
48namespace internal
49{
62 template <int dim>
63 constexpr ReferenceCell<dim>
64 make_reference_cell_from_int(const std::uint8_t kind);
65} // namespace internal
66
67
68
75enum class IsotropicRefinementChoice : std::uint8_t
76{
80 isotropic_refinement = 0,
84 cut_tet_68 = 1,
88 cut_tet_57 = 2,
92 cut_tet_49 = 3,
93};
94
95
96
127template <int dim>
128class ReferenceCell
129{
130public:
140 DEAL_II_DEPRECATED_WITH_COMMENT(
141 "Use ReferenceCells::n_vertices_to_reference_cell() instead")
142 static ReferenceCell
143 n_vertices_to_type(const int dim_, const unsigned int n_vertices);
144
153 constexpr ReferenceCell();
154
165 constexpr bool
166 is_hyper_cube() const;
167
171 constexpr bool
172 is_simplex() const;
173
179 constexpr unsigned int
180 get_dimension() const;
181
201 double
202 d_linear_shape_function(const Point<dim> &xi, const unsigned int i) const;
203
208 Tensor<1, dim>
209 d_linear_shape_function_gradient(const Point<dim> &xi,
210 const unsigned int i) const;
211
221 template <int spacedim = dim>
222 std::unique_ptr<Mapping<dim, spacedim>>
223 get_default_mapping(const unsigned int degree) const;
224
236 template <int spacedim = dim>
237 const Mapping<dim, spacedim> &
238 get_default_linear_mapping() const;
239
248 Quadrature<dim>
249 get_gauss_type_quadrature(const unsigned n_points_1d) const;
250
260 Quadrature<dim>
261 get_midpoint_quadrature() const;
262
276 const Quadrature<dim> &
277 get_nodal_type_quadrature() const;
278
293 constexpr unsigned int
294 n_vertices() const;
295
300 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
301 vertex_indices() const;
302
310 Point<dim>
311 vertex(const unsigned int v) const;
312
318 constexpr unsigned int
319 n_lines() const;
320
325 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
326 line_indices() const;
327
333 constexpr unsigned int
334 n_faces() const;
335
340 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
341 face_indices() const;
342
348 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
349 face_indices_by_type(const ReferenceCell<dim - 1> &face_ref_type) const;
350
354 constexpr unsigned int
355 n_isotropic_refinement_choices() const;
356
360 IsotropicRefinementChoice
361 get_isotropic_refinement_choice(const unsigned int ref_choice) const;
362
374 constexpr unsigned int
375 n_isotropic_children() const;
376
381 unsigned int
382 n_children(const RefinementCase<dim> ref_case =
383 RefinementCase<dim>::isotropic_refinement) const;
384
389 std_cxx20::ranges::iota_view<unsigned int, unsigned int>
390 isotropic_child_indices() const;
391
411 std::pair<unsigned int, RefinementCase<dim - 1>>
412 equivalent_refinement_case(
413 const types::geometric_orientation combined_face_orientation,
414 const internal::SubfaceCase<dim> subface_case,
415 const unsigned int subface_no) const;
416
424 ArrayView<const RefinementCase<dim>>
425 refinement_cases() const;
426
439 constexpr ReferenceCell<dim - 1>
440 face_reference_cell(const unsigned int face_index) const;
441
491 unsigned int
492 child_cell_on_face(
493 const unsigned int face,
494 const unsigned int subface,
495 const types::geometric_orientation combined_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
551 unsigned int
552 face_to_cell_lines(const unsigned int face,
553 const unsigned int line,
554 const types::geometric_orientation face_orientation) const;
555
572 unsigned int
573 face_to_cell_vertices(
574 const unsigned int face,
575 const unsigned int vertex,
576 const types::geometric_orientation face_orientation) const;
577
598 unsigned int
599 face_and_line_to_cell_vertices(
600 const unsigned int face,
601 const unsigned int line,
602 const unsigned int vertex,
603 const types::geometric_orientation face_orientation) const;
604
624 Point<dim>
625 face_vertex_location(const unsigned int face,
626 const unsigned int vertex) const;
627
641 Point<dim>
642 subface_vertex_location(
643 const unsigned int face_no,
644 const unsigned int subface_no,
645 const unsigned int subface_vertex_no,
646 const RefinementCase<dim - 1> face_refinement_case) const;
647
652 unsigned int
653 standard_to_real_face_vertex(
654 const unsigned int vertex,
655 const unsigned int face,
656 const types::geometric_orientation face_orientation) const;
657
662 unsigned int
663 standard_to_real_face_line(
664 const unsigned int line,
665 const unsigned int face,
666 const types::geometric_orientation face_orientation) const;
667
690 types::geometric_orientation
691 face_to_cell_line_orientation(
692 const unsigned int face_line_no,
693 const unsigned int face_no,
694 const types::geometric_orientation face_orientation,
695 const types::geometric_orientation line_orientation) const;
696
722 unsigned int
723 opposite_face_index(const unsigned int face_no) const;
724
748 double
749 volume() const;
750
766 double
767 face_measure(const unsigned int face_no) const;
768
781 Point<dim>
782 barycenter() const;
783
803 bool
804 contains_point(const Point<dim> &p, const double tolerance = 0) const;
805
810 Point<dim>
811 closest_point(const Point<dim> &p) const;
812
820 Tensor<1, dim>
821 face_tangent_vector(const unsigned int face_no, const unsigned int i) const;
822
826 Tensor<1, dim>
827 face_normal_vector(const unsigned int face_no) const;
828
834 unsigned int
835 n_face_orientations(const unsigned int face_no) const;
836
853 template <typename T, std::size_t N>
854 DEAL_II_DEPRECATED types::geometric_orientation
855 compute_orientation(const std::array<T, N> &vertices_0,
856 const std::array<T, N> &vertices_1) const;
857
882 template <typename T>
883 types::geometric_orientation
884 get_combined_orientation(const ArrayView<const T> &vertices_0,
885 const ArrayView<const T> &vertices_1) const;
886
898 template <typename T>
899 boost::container::small_vector<T, 8>
900 permute_by_combined_orientation(
901 const ArrayView<const T> &vertices,
902 const types::geometric_orientation orientation) const;
903
909 types::geometric_orientation
910 get_inverse_combined_orientation(
911 const types::geometric_orientation orientation) const;
912
916 ArrayView<const unsigned int>
917 faces_for_given_vertex(const unsigned int vertex_index) const;
918
930 const ndarray<unsigned int, 13, 4> &
931 new_isotropic_child_face_lines(const unsigned int refinement_choice) const;
932
946 const ndarray<unsigned int, 13, 4, 2> &
947 new_isotropic_child_face_line_vertices(
948 const unsigned int refinement_choice) const;
949
961 const ndarray<unsigned int, 10, 6> &
962 new_isotropic_child_cell_faces(const unsigned int refinement_choice) const;
963
975 const ndarray<unsigned int, 10, 8> &
976 new_isotropic_child_cell_vertices(const unsigned int refinement_choice) const;
977
990 unsigned int
991 exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const;
992
996 unsigned int
997 exodusii_face_to_deal_face(const unsigned int face_n) const;
998
1002 unsigned int
1003 ucd_vertex_to_deal_vertex(const unsigned int vertex_n) const;
1004
1008 unsigned int
1009 unv_vertex_to_deal_vertex(const unsigned int vertex_n) const;
1010
1014 unsigned int
1015 vtk_linear_type() const;
1016
1021 unsigned int
1022 vtk_quadratic_type() const;
1023
1028 unsigned int
1029 vtk_lagrange_type() const;
1030
1043 unsigned int
1044 vtk_lexicographic_to_node_index(
1045 const std::array<unsigned, dim> &node_indices,
1046 const std::array<unsigned, dim> &nodes_per_direction,
1047 const bool legacy_format) const;
1048
1052 unsigned int
1053 vtk_vertex_to_deal_vertex(const unsigned int vertex_index) const;
1054
1058 unsigned int
1059 gmsh_element_type() const;
1060
1074 std::string
1075 to_string() const;
1076
1080 constexpr operator std::uint8_t() const;
1081
1085 constexpr bool
1086 operator==(const ReferenceCell<dim> &type) const;
1087
1095 template <int dim2>
1096 constexpr bool
1097 operator==(const ReferenceCell<dim2> &type) const;
1098
1102 constexpr bool
1103 operator!=(const ReferenceCell<dim> &type) const;
1104
1112 template <int dim2>
1113 constexpr bool
1114 operator!=(const ReferenceCell<dim2> &type) const;
1115
1121 template <class Archive>
1122 void
1123 serialize(Archive &archive, const unsigned int /*version*/);
1124
1128 static constexpr std::size_t
1129 memory_consumption();
1130
1135private:
1139 std::uint8_t kind;
1140
1147 explicit constexpr ReferenceCell(const std::uint8_t kind);
1148
1154 static constexpr ndarray<unsigned int, 1, 1> vertex_vertex_permutations = {
1155 {{{0}}}};
1156
1163
1174
1179 static constexpr ndarray<unsigned int, 8, 4>
1180 quadrilateral_vertex_permutations = {{
1181 {{0, 1, 2, 3}},
1182 {{0, 2, 1, 3}},
1183 {{2, 0, 3, 1}},
1184 {{2, 3, 0, 1}},
1185 {{3, 2, 1, 0}},
1186 {{3, 1, 2, 0}},
1187 {{1, 3, 0, 2}},
1188 {{1, 0, 3, 2}},
1189 }};
1190
1194 static constexpr ndarray<unsigned int, 6, 3> triangle_line_permutations = {
1195 {{{0, 1, 2}},
1196 {{2, 1, 0}},
1197 {{2, 0, 1}},
1198 {{1, 0, 2}},
1199 {{1, 2, 0}},
1200 {{0, 2, 1}}}};
1201
1206 static constexpr ndarray<unsigned int, 8, 4> quadrilateral_line_permutations =
1207 {{
1208 {{0, 1, 2, 3}},
1209 {{2, 3, 0, 1}},
1210 {{3, 2, 0, 1}},
1211 {{0, 1, 3, 2}},
1212 {{1, 0, 3, 2}},
1213 {{3, 2, 1, 0}},
1214 {{2, 3, 1, 0}},
1215 {{1, 0, 2, 3}},
1216 }};
1217
1222 template <int dim_>
1223 friend constexpr ReferenceCell<dim_>
1224 internal::make_reference_cell_from_int(const std::uint8_t);
1225
1226 template <int dim_>
1227 friend std::ostream &
1228 operator<<(std::ostream &out, const ReferenceCell<dim_> &reference_cell);
1229
1230 template <int dim_>
1231 friend std::istream &
1232 operator>>(std::istream &in, ReferenceCell<dim_> &reference_cell);
1233};
1234
1235
1244template <int dim>
1245std::ostream &
1246operator<<(std::ostream &out, const ReferenceCell<dim> &reference_cell);
1247
1255template <int dim>
1256std::istream &
1257operator>>(std::istream &in, ReferenceCell<dim> &reference_cell);
1258
1259
1260
1261template <int dim>
1262inline constexpr ReferenceCell<dim>::ReferenceCell(const std::uint8_t kind)
1263 : kind(kind)
1264{}
1265
1266
1267
1268template <int dim>
1269inline constexpr ReferenceCell<dim>::operator std::uint8_t() const
1270{
1271 return kind;
1272}
1273
1274
1275
1276template <int dim>
1277inline constexpr bool
1278ReferenceCell<dim>::operator==(const ReferenceCell<dim> &type) const
1279{
1280 return kind == type.kind;
1281}
1282
1283
1284
1285template <int dim>
1286template <int dim2>
1287inline constexpr bool
1288ReferenceCell<dim>::operator==(const ReferenceCell<dim2> &) const
1289{
1290 return false;
1291}
1292
1293
1294
1295template <int dim>
1296inline constexpr bool
1297ReferenceCell<dim>::operator!=(const ReferenceCell<dim> &type) const
1298{
1299 return kind != type.kind;
1300}
1301
1302
1303
1304template <int dim>
1305template <int dim2>
1306inline constexpr bool
1307ReferenceCell<dim>::operator!=(const ReferenceCell<dim2> &) const
1308{
1309 return true;
1310}
1311
1312
1313
1314namespace internal
1315{
1316 template <int dim>
1317 inline constexpr ReferenceCell<dim>
1318 make_reference_cell_from_int(const std::uint8_t kind)
1319 {
1320#ifndef DEAL_II_CXX14_CONSTEXPR_BUG
1321 // Make sure these are the only indices from which objects can be
1322 // created.
1323 Assert((kind == std::numeric_limits<std::uint8_t>::max()) || (kind < 8),
1324 ExcInternalError());
1325#endif
1326
1327 // Call the private constructor, which we can from here because this
1328 // function is a 'friend'.
1329 return ReferenceCell<dim>(kind);
1330 }
1331} // namespace internal
1332
1333
1334
1342namespace ReferenceCells
1343{
1344 constexpr ReferenceCell<0> Vertex =
1345 internal::make_reference_cell_from_int<0>(0);
1346 constexpr ReferenceCell<1> Line =
1347 internal::make_reference_cell_from_int<1>(1);
1348 constexpr ReferenceCell<2> Triangle =
1349 internal::make_reference_cell_from_int<2>(2);
1350 constexpr ReferenceCell<2> Quadrilateral =
1351 internal::make_reference_cell_from_int<2>(3);
1352 constexpr ReferenceCell<3> Tetrahedron =
1353 internal::make_reference_cell_from_int<3>(4);
1354 constexpr ReferenceCell<3> Pyramid =
1355 internal::make_reference_cell_from_int<3>(5);
1356 constexpr ReferenceCell<3> Wedge =
1357 internal::make_reference_cell_from_int<3>(6);
1358 constexpr ReferenceCell<3> Hexahedron =
1359 internal::make_reference_cell_from_int<3>(7);
1360
1361 template <int dim>
1362 constexpr ReferenceCell<dim> Invalid =
1363 internal::make_reference_cell_from_int<dim>(
1364 std::numeric_limits<std::uint8_t>::max());
1365
1370 template <int dim>
1371 constexpr const std::vector<ReferenceCell<dim>>
1372 get_reference_cells_in_dim();
1373
1379 template <int dim>
1380 constexpr const ReferenceCell<dim> &
1381 get_simplex();
1382
1388 template <int dim>
1389 constexpr const ReferenceCell<dim> &
1390 get_hypercube();
1391
1396 template <int dim>
1397 constexpr unsigned int
1398 n_reference_cells_in_dim();
1399
1412 template <int dim>
1413 ReferenceCell<dim>
1414 n_vertices_to_reference_cell(const unsigned int n_vertices);
1415
1423 template <int dim>
1424 constexpr unsigned int
1425 max_n_vertices();
1426
1434 template <int dim>
1435 constexpr unsigned int
1436 max_n_lines();
1437
1453 template <int dim>
1454 constexpr unsigned int
1455 max_n_faces();
1456
1466 template <int dim>
1467 constexpr unsigned int
1468 max_n_children();
1469} // namespace ReferenceCells
1470
1471
1472
1473template <int dim>
1474inline constexpr ReferenceCell<dim>::ReferenceCell()
1475 : ReferenceCell(ReferenceCells::Invalid<dim>)
1476{}
1477
1478
1479
1480template <int dim>
1481template <class Archive>
1482inline void
1483ReferenceCell<dim>::serialize(Archive &archive, const unsigned int /*version*/)
1484{
1485 archive &kind;
1486}
1487
1488
1489
1490template <int dim>
1491inline constexpr std::size_t
1492ReferenceCell<dim>::memory_consumption()
1493{
1494 return sizeof(ReferenceCells::Invalid<dim>);
1495}
1496
1497
1498
1499template <int dim>
1500inline ArrayView<const unsigned int>
1501ReferenceCell<dim>::faces_for_given_vertex(const unsigned int vertex) const
1502{
1503 AssertIndexRange(vertex, n_vertices());
1504 if constexpr (dim == 0)
1505 DEAL_II_NOT_IMPLEMENTED();
1506 else if constexpr (dim == 1)
1507 return {&GeometryInfo<1>::vertex_to_face[vertex][0], 1};
1508 else if constexpr (dim == 2)
1509 {
1510 switch (this->kind)
1511 {
1512 case ReferenceCells::Quadrilateral:
1513 return {&GeometryInfo<2>::vertex_to_face[vertex][0], 2};
1514 case ReferenceCells::Triangle:
1515 {
1516 static constexpr ndarray<unsigned int, 3, 2> table = {
1517 {{{0, 2}}, {{0, 1}}, {{1, 2}}}};
1518 return table[vertex];
1519 }
1520 default:
1521 DEAL_II_NOT_IMPLEMENTED();
1522 }
1523 }
1524 else if constexpr (dim == 3)
1525 {
1526 switch (this->kind)
1527 {
1528 case ReferenceCells::Tetrahedron:
1529 {
1530 static constexpr ndarray<unsigned int, 4, 3> table = {
1531 {{{0, 1, 2}}, {{0, 1, 3}}, {{0, 2, 3}}, {{1, 2, 3}}}};
1532
1533 return table[vertex];
1534 }
1535 case ReferenceCells::Pyramid:
1536 {
1537 static constexpr unsigned int X = numbers::invalid_unsigned_int;
1538 static constexpr ndarray<unsigned int, 5, 4> table = {
1539 {{{0, 1, 3, X}},
1540 {{0, 2, 3, X}},
1541 {{0, 1, 4, X}},
1542 {{0, 2, 4, X}},
1543 {{1, 2, 3, 4}}}};
1544
1545 return {&table[vertex][0], vertex == 4 ? 4u : 3u};
1546 }
1547 case ReferenceCells::Wedge:
1548 {
1549 AssertIndexRange(vertex, 6);
1550 static constexpr ndarray<unsigned int, 6, 3> table = {
1551 {{{0, 2, 4}},
1552 {{0, 2, 3}},
1553 {{0, 3, 4}},
1554 {{1, 2, 4}},
1555 {{1, 2, 3}},
1556 {{1, 3, 4}}}};
1557
1558 return table[vertex];
1559 }
1560 case ReferenceCells::Hexahedron:
1561 return {&GeometryInfo<3>::vertex_to_face[vertex][0], 3};
1562 default:
1563 DEAL_II_NOT_IMPLEMENTED();
1564 }
1565 }
1566 else
1567 DEAL_II_NOT_IMPLEMENTED();
1568
1569 return {};
1570}
1571
1572
1573
1574template <int dim>
1575const ndarray<unsigned int, 13, 4> &
1576ReferenceCell<dim>::new_isotropic_child_face_lines(
1577 const unsigned int refinement_choice) const
1578{
1579 AssertIndexRange(refinement_choice, n_isotropic_refinement_choices());
1580
1581 const unsigned int X = numbers::invalid_unsigned_int;
1582 if constexpr (dim == 3)
1583 {
1584 switch (this->kind)
1585 {
1586 case ReferenceCells::Tetrahedron:
1587 switch (refinement_choice)
1588 {
1589 case 0:
1590 {
1591 // new line is (6,8)
1592 static constexpr ndarray<unsigned int, 13, 4> table_68 = {
1593 {{{2, 3, 8, X}},
1594 {{0, 9, 5, X}},
1595 {{1, 6, 11, X}},
1596 {{4, 10, 7, X}},
1597 {{2, 12, 5, X}},
1598 {{1, 9, 12, X}},
1599 {{4, 8, 12, X}},
1600 {{6, 12, 10, X}},
1601 {{X, X, X, X}},
1602 {{X, X, X, X}},
1603 {{X, X, X, X}},
1604 {{X, X, X, X}},
1605 {{X, X, X, X}}}};
1606 return table_68;
1607 }
1608 case 1:
1609 {
1610 // new line is (5,7)
1611 static constexpr ndarray<unsigned int, 13, 4> table_57 = {
1612 {{{2, 3, 8, X}},
1613 {{0, 9, 5, X}},
1614 {{1, 6, 11, X}},
1615 {{4, 10, 7, X}},
1616 {{0, 3, 12, X}},
1617 {{1, 12, 8, X}},
1618 {{4, 12, 9, X}},
1619 {{7, 11, 12, X}},
1620 {{X, X, X, X}},
1621 {{X, X, X, X}},
1622 {{X, X, X, X}},
1623 {{X, X, X, X}},
1624 {{X, X, X, X}}}};
1625 return table_57;
1626 }
1627 case 2:
1628 {
1629 // new line is (4,9)
1630 static constexpr ndarray<unsigned int, 13, 4> table_49 = {
1631 {{{2, 3, 8, X}},
1632 {{0, 9, 5, X}},
1633 {{1, 6, 11, X}},
1634 {{4, 10, 7, X}},
1635 {{0, 12, 11, X}},
1636 {{2, 6, 12, X}},
1637 {{3, 12, 7, X}},
1638 {{5, 10, 12, X}},
1639 {{X, X, X, X}},
1640 {{X, X, X, X}},
1641 {{X, X, X, X}},
1642 {{X, X, X, X}},
1643 {{X, X, X, X}}}};
1644 return table_49;
1645 }
1646 default:
1647 DEAL_II_ASSERT_UNREACHABLE();
1648 }
1649
1650 case ReferenceCells::Pyramid:
1651 {
1652 static constexpr ndarray<unsigned int, 13, 4>
1653 new_quad_lines_pyramid = {
1654 {{{0, 10, 16, X}}, // child 0
1655 {{2, 16, 6, X}},
1656 {{0, 17, 12, X}}, // child 1
1657 {{3, 7, 17, X}},
1658 {{1, 18, 15, X}}, // child 2
1659 {{2, 4, 18, X}},
1660 {{1, 13, 19, X}}, // child 3
1661 {{3, 19, 9, X}},
1662 {{5, 8, 11, 14}}, // child 9 (top)
1663 {{11, 17, 16, X}}, // child 8 (upside down)
1664 {{14, 18, 19, X}},
1665 {{5, 16, 18, X}},
1666 {{8, 19, 17, X}}}};
1667 return new_quad_lines_pyramid;
1668 }
1669
1670 case ReferenceCells::Wedge:
1671 {
1672 static constexpr ndarray<unsigned int, 13, 4>
1673 new_quad_lines_wedge = {{{{2, 20, 11, X}}, // mid tri
1674 {{3, 6, 18, X}},
1675 {{19, 7, 10, X}},
1676 {{18, 19, 20, X}},
1677 {{4, 0, 14, 18}}, // lower
1678 {{0, 8, 12, 20}},
1679 {{8, 4, 13, 19}},
1680 {{5, 1, 18, 15}}, // upper
1681 {{1, 9, 20, 17}},
1682 {{9, 5, 19, 16}},
1683 {{X, X, X, X}},
1684 {{X, X, X, X}},
1685 {{X, X, X, X}}}};
1686 return new_quad_lines_wedge;
1687 }
1688
1689 case ReferenceCells::Hexahedron:
1690 {
1691 static constexpr ndarray<unsigned int, 13, 4> new_quad_lines_hex =
1692 {{{{10, 28, 16, 24}},
1693 {{28, 14, 17, 25}},
1694 {{11, 29, 24, 20}},
1695 {{29, 15, 25, 21}},
1696 {{18, 26, 0, 28}},
1697 {{26, 22, 1, 29}},
1698 {{19, 27, 28, 4}},
1699 {{27, 23, 29, 5}},
1700 {{2, 24, 8, 26}},
1701 {{24, 6, 9, 27}},
1702 {{3, 25, 26, 12}},
1703 {{25, 7, 27, 13}},
1704 {{X, X, X, X}}}};
1705 return new_quad_lines_hex;
1706 }
1707
1708 default:
1709 DEAL_II_ASSERT_UNREACHABLE();
1710 }
1711 }
1712 else
1713 // We should never get here except in the 3d case.
1714 DEAL_II_ASSERT_UNREACHABLE();
1715
1716 static constexpr ndarray<unsigned int, 13, 4> empty{};
1717 return empty;
1718}
1719
1720
1721
1722template <int dim>
1723const ndarray<unsigned int, 13, 4, 2> &
1724ReferenceCell<dim>::new_isotropic_child_face_line_vertices(
1725 const unsigned int refinement_choice) const
1726{
1727 AssertIndexRange(refinement_choice, n_isotropic_refinement_choices());
1728
1729 const unsigned int X = numbers::invalid_unsigned_int;
1730 if constexpr (dim == 3)
1731 {
1732 switch (this->kind)
1733 {
1734 case ReferenceCells::Tetrahedron:
1735 switch (refinement_choice)
1736 {
1737 case 0:
1738 {
1739 // new line is (6, 8)
1740 static constexpr ndarray<unsigned int, 13, 4, 2> table_68 =
1741 {{{{{{6, 4}}, {{4, 7}}, {{7, 6}}, {{X, X}}}},
1742 {{{{4, 5}}, {{5, 8}}, {{8, 4}}, {{X, X}}}},
1743 {{{{5, 6}}, {{6, 9}}, {{9, 5}}, {{X, X}}}},
1744 {{{{7, 8}}, {{8, 9}}, {{9, 7}}, {{X, X}}}},
1745 {{{{4, 6}}, {{6, 8}}, {{8, 4}}, {{X, X}}}},
1746 {{{{6, 5}}, {{5, 8}}, {{8, 6}}, {{X, X}}}},
1747 {{{{8, 7}}, {{7, 6}}, {{6, 8}}, {{X, X}}}},
1748 {{{{9, 6}}, {{6, 8}}, {{8, 9}}, {{X, X}}}},
1749 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1750 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1751 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1752 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1753 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}}}};
1754 return table_68;
1755 }
1756 case 1:
1757 {
1758 // new line is (5, 7)
1759 static constexpr ndarray<unsigned int, 13, 4, 2> table_57 =
1760 {{{{{{6, 4}}, {{4, 7}}, {{7, 6}}, {{X, X}}}},
1761 {{{{4, 5}}, {{5, 8}}, {{8, 4}}, {{X, X}}}},
1762 {{{{5, 6}}, {{6, 9}}, {{9, 5}}, {{X, X}}}},
1763 {{{{7, 8}}, {{8, 9}}, {{9, 7}}, {{X, X}}}},
1764 {{{{5, 4}}, {{4, 7}}, {{7, 5}}, {{X, X}}}},
1765 {{{{6, 5}}, {{5, 7}}, {{7, 6}}, {{X, X}}}},
1766 {{{{8, 7}}, {{7, 5}}, {{5, 8}}, {{X, X}}}},
1767 {{{{7, 9}}, {{9, 5}}, {{5, 7}}, {{X, X}}}},
1768 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1769 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1770 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1771 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1772 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}}}};
1773 return table_57;
1774 }
1775 case 2:
1776 {
1777 // new line is (4, 9)
1778 static constexpr ndarray<unsigned int, 13, 4, 2> table_49 =
1779 {{{{{{6, 4}}, {{4, 7}}, {{7, 6}}, {{X, X}}}},
1780 {{{{4, 5}}, {{5, 8}}, {{8, 4}}, {{X, X}}}},
1781 {{{{5, 6}}, {{6, 9}}, {{9, 5}}, {{X, X}}}},
1782 {{{{7, 8}}, {{8, 9}}, {{9, 7}}, {{X, X}}}},
1783 {{{{5, 4}}, {{4, 9}}, {{9, 5}}, {{X, X}}}},
1784 {{{{4, 6}}, {{6, 9}}, {{9, 4}}, {{X, X}}}},
1785 {{{{7, 4}}, {{4, 9}}, {{9, 7}}, {{X, X}}}},
1786 {{{{4, 8}}, {{8, 9}}, {{9, 4}}, {{X, X}}}},
1787 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1788 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1789 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1790 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1791 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}}}};
1792 return table_49;
1793 }
1794 default:
1795 DEAL_II_ASSERT_UNREACHABLE();
1796 }
1797
1798 case ReferenceCells::Pyramid:
1799 {
1800 static constexpr ndarray<unsigned int, 13, 4, 2>
1801 quad_lines_vertices_pyramid = {
1802 {// child 0
1803 {{{{13, 7}}, {{7, 9}}, {{9, 13}}, {{X, X}}}},
1804 {{{{5, 13}}, {{13, 9}}, {{9, 5}}, {{X, X}}}},
1805 // child 1
1806 {{{{7, 13}}, {{13, 10}}, {{10, 7}}, {{X, X}}}},
1807 {{{{13, 6}}, {{6, 10}}, {{10, 13}}, {{X, X}}}},
1808 // child 2
1809 {{{{8, 13}}, {{13, 11}}, {{11, 8}}, {{X, X}}}},
1810 {{{{13, 5}}, {{5, 11}}, {{11, 13}}, {{X, X}}}},
1811 // child 3
1812 {{{{13, 8}}, {{8, 12}}, {{12, 13}}, {{X, X}}}},
1813 {{{{6, 13}}, {{13, 12}}, {{12, 6}}, {{X, X}}}},
1814 // child 9 (top)
1815 {{{{9, 11}}, {{10, 12}}, {{9, 10}}, {{11, 12}}}},
1816 // child 8 (upside down)
1817 {{{{9, 10}}, {{10, 13}}, {{13, 9}}, {{X, X}}}},
1818 {{{{12, 11}}, {{11, 13}}, {{13, 12}}, {{X, X}}}},
1819 {{{{11, 9}}, {{9, 13}}, {{13, 11}}, {{X, X}}}},
1820 {{{{10, 12}}, {{12, 13}}, {{13, 10}}, {{X, X}}}}}};
1821 return quad_lines_vertices_pyramid;
1822 }
1823
1824 case ReferenceCells::Wedge:
1825 {
1826 static constexpr ndarray<unsigned int, 13, 4, 2>
1827 quad_lines_vertices_wedge = {
1828 {// mid tri
1829 {{{{12, 15}}, {{15, 17}}, {{17, 12}}, {{X, X}}}},
1830 {{{{15, 13}}, {{13, 16}}, {{16, 15}}, {{X, X}}}},
1831 {{{{17, 16}}, {{16, 14}}, {{14, 17}}, {{X, X}}}},
1832 {{{{15, 16}}, {{16, 17}}, {{17, 15}}, {{X, X}}}},
1833 // lower
1834 {{{{7, 16}}, {{6, 15}}, {{7, 6}}, {{16, 15}}}},
1835 {{{{6, 15}}, {{8, 17}}, {{6, 8}}, {{15, 17}}}},
1836 {{{{8, 17}}, {{7, 16}}, {{8, 7}}, {{17, 16}}}},
1837 // upper
1838 {{{{16, 10}}, {{15, 9}}, {{16, 15}}, {{10, 9}}}},
1839 {{{{15, 9}}, {{17, 11}}, {{15, 17}}, {{9, 11}}}},
1840 {{{{17, 11}}, {{16, 10}}, {{17, 16}}, {{11, 10}}}},
1841 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1842 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}},
1843 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}}}};
1844 return quad_lines_vertices_wedge;
1845 }
1846
1847 case ReferenceCells::Hexahedron:
1848 {
1849 static constexpr ndarray<unsigned int, 13, 4, 2>
1850 quad_lines_vertices_hex = {
1851 {{{{{10, 22}}, {{24, 26}}, {{10, 24}}, {{22, 26}}}},
1852 {{{{24, 26}}, {{11, 23}}, {{24, 11}}, {{26, 23}}}},
1853 {{{{22, 14}}, {{26, 25}}, {{22, 26}}, {{14, 25}}}},
1854 {{{{26, 25}}, {{23, 15}}, {{26, 23}}, {{25, 15}}}},
1855 {{{{8, 24}}, {{20, 26}}, {{8, 20}}, {{24, 26}}}},
1856 {{{{20, 26}}, {{12, 25}}, {{20, 12}}, {{26, 25}}}},
1857 {{{{24, 9}}, {{26, 21}}, {{24, 26}}, {{9, 21}}}},
1858 {{{{26, 21}}, {{25, 13}}, {{26, 25}}, {{21, 13}}}},
1859 {{{{16, 20}}, {{22, 26}}, {{16, 22}}, {{20, 26}}}},
1860 {{{{22, 26}}, {{17, 21}}, {{22, 17}}, {{26, 21}}}},
1861 {{{{20, 18}}, {{26, 23}}, {{20, 26}}, {{18, 23}}}},
1862 {{{{26, 23}}, {{21, 19}}, {{26, 21}}, {{23, 19}}}},
1863 {{{{X, X}}, {{X, X}}, {{X, X}}, {{X, X}}}}}};
1864 return quad_lines_vertices_hex;
1865 }
1866 default:
1867 DEAL_II_ASSERT_UNREACHABLE();
1868 }
1869 }
1870 else
1871 // We should never get here except in the 3d case.
1872 DEAL_II_ASSERT_UNREACHABLE();
1873
1874 static constexpr ndarray<unsigned int, 13, 4, 2> empty{};
1875 return empty;
1876}
1877
1878
1879
1880template <int dim>
1881const ndarray<unsigned int, 10, 6> &
1882ReferenceCell<dim>::new_isotropic_child_cell_faces(
1883 const unsigned int refinement_choice) const
1884{
1885 AssertIndexRange(refinement_choice, n_isotropic_refinement_choices());
1886
1887 const unsigned int X = numbers::invalid_unsigned_int;
1888 if constexpr (dim == 3)
1889 {
1890 switch (this->kind)
1891 {
1892 case ReferenceCells::Tetrahedron:
1893 switch (refinement_choice)
1894 {
1895 case 0:
1896 {
1897 // new line is (6, 8)
1898 static constexpr ndarray<unsigned int, 10, 6> table_68 = {{
1899 {{8, 13, 16, 0, X, X}},
1900 {{9, 12, 1, 21, X, X}},
1901 {{10, 2, 17, 20, X, X}},
1902 {{3, 14, 18, 22, X, X}},
1903 {{11, 1, 4, 5, X, X}},
1904 {{15, 0, 4, 6, X, X}},
1905 {{19, 7, 6, 3, X, X}},
1906 {{23, 5, 2, 7, X, X}},
1907 {{X, X, X, X, X, X}},
1908 {{X, X, X, X, X, X}},
1909 }};
1910 return table_68;
1911 }
1912 case 1:
1913 {
1914 // new line is (5, 7)
1915 static constexpr ndarray<unsigned int, 10, 6> table_57 = {{
1916 {{8, 13, 16, 0, X, X}},
1917 {{9, 12, 1, 21, X, X}},
1918 {{10, 2, 17, 20, X, X}},
1919 {{3, 14, 18, 22, X, X}},
1920 {{11, 4, 0, 5, X, X}},
1921 {{15, 4, 1, 6, X, X}},
1922 {{19, 2, 5, 7, X, X}},
1923 {{23, 6, 7, 3, X, X}},
1924 {{X, X, X, X, X, X}},
1925 {{X, X, X, X, X, X}},
1926 }};
1927 return table_57;
1928 }
1929 case 2:
1930 {
1931 // new line is (4, 9)
1932 static constexpr ndarray<unsigned int, 10, 6> table_49 = {{
1933 {{8, 13, 16, 0, X, X}},
1934 {{9, 12, 1, 21, X, X}},
1935 {{10, 2, 17, 20, X, X}},
1936 {{3, 14, 18, 22, X, X}},
1937 {{11, 4, 5, 2, X, X}},
1938 {{15, 6, 7, 3, X, X}},
1939 {{19, 5, 0, 6, X, X}},
1940 {{23, 1, 4, 7, X, X}},
1941 {{X, X, X, X, X, X}},
1942 {{X, X, X, X, X, X}},
1943 }};
1944 return table_49;
1945 }
1946 default:
1947 DEAL_II_ASSERT_UNREACHABLE();
1948 }
1949
1950 case ReferenceCells::Pyramid:
1951 {
1952 static constexpr ndarray<unsigned int, 10, 6> cell_quads_pyramid =
1953 {{
1954 {{13, 17, 0, 26, 1, X}}, // bottom pyramids
1955 {{14, 2, 22, 25, 3, X}}, //
1956 {{15, 18, 4, 5, 29, X}}, //
1957 {{16, 6, 21, 7, 30, X}}, //
1958 {{5, 1, 20, 11, X, X}}, // bottom wedges
1959 {{7, 3, 12, 24, X, X}}, //
1960 {{0, 2, 28, 9, X, X}}, //
1961 {{4, 6, 10, 32, X, X}}, //
1962 {{8, 9, 10, 11, 12, X}}, // upside down pyramid
1963 {{8, 19, 23, 27, 31, X}} // top pyramid
1964 }};
1965 return cell_quads_pyramid;
1966 }
1967
1968 case ReferenceCells::Wedge:
1969 {
1970 // Considerations that went into creating this table
1971 // - Keep aspect ratio of 1.0 if parent had 1.0 (this is why the
1972 // centre wedge is rotated by 180 degrees compared to the
1973 // parent)
1974 // - Order of children is similar to the order of a reference
1975 // triangle in the x-y-plane (i.e. like the top triangle of the
1976 // parent)
1977 static constexpr ndarray<unsigned int, 10, 6> cell_quads_wedge = {
1978 {
1979 {{23, 0, 10, 5, 19, X}}, // bottom children
1980 {{22, 1, 11, 14, 4, X}}, //
1981 {{24, 2, 6, 15, 18, X}}, //
1982 {{25, 3, 6, 5, 4, X}}, //
1983 {{0, 26, 12, 8, 21, X}}, // top children
1984 {{1, 27, 13, 16, 7, X}}, //
1985 {{2, 28, 9, 17, 20, X}}, //
1986 {{3, 29, 9, 8, 7, X}}, //
1987 {{X, X, X, X, X, X}},
1988 {{X, X, X, X, X, X}},
1989 }};
1990 return cell_quads_wedge;
1991 }
1992
1993 case ReferenceCells::Hexahedron:
1994 {
1995 static constexpr ndarray<unsigned int, 10, 6> cell_quads_hex = {{
1996 {{12, 0, 20, 4, 28, 8}}, // bottom children
1997 {{0, 16, 22, 6, 29, 9}}, //
1998 {{13, 1, 4, 24, 30, 10}}, //
1999 {{1, 17, 6, 26, 31, 11}}, //
2000 {{14, 2, 21, 5, 8, 32}}, // top children
2001 {{2, 18, 23, 7, 9, 33}}, //
2002 {{15, 3, 5, 25, 10, 34}}, //
2003 {{3, 19, 7, 27, 11, 35}}, //
2004 {{X, X, X, X, X, X}},
2005 {{X, X, X, X, X, X}},
2006 }};
2007 return cell_quads_hex;
2008 }
2009 default:
2010 DEAL_II_ASSERT_UNREACHABLE();
2011 }
2012 }
2013 else
2014 // We should never get here except in the 3d case.
2015 DEAL_II_ASSERT_UNREACHABLE();
2016
2017 static constexpr ndarray<unsigned int, 10, 6> empty{};
2018 return empty;
2019}
2020
2021
2022
2023template <int dim>
2024const ndarray<unsigned int, 10, 8> &
2025ReferenceCell<dim>::new_isotropic_child_cell_vertices(
2026 const unsigned int refinement_choice) const
2027{
2028 AssertIndexRange(refinement_choice, n_isotropic_refinement_choices());
2029
2030 constexpr unsigned int X = numbers::invalid_unsigned_int;
2031 if constexpr (dim == 3)
2032 {
2033 switch (this->kind)
2034 {
2035 case ReferenceCells::Tetrahedron:
2036 switch (refinement_choice)
2037 {
2038 case 0:
2039 {
2040 // new line is (6,8)
2041 static constexpr ndarray<unsigned int, 10, 8> table_68 = {{
2042 {{0, 4, 6, 7, X, X, X, X}},
2043 {{4, 1, 5, 8, X, X, X, X}},
2044 {{6, 5, 2, 9, X, X, X, X}},
2045 {{7, 8, 9, 3, X, X, X, X}},
2046 {{4, 5, 6, 8, X, X, X, X}},
2047 {{4, 7, 8, 6, X, X, X, X}},
2048 {{6, 9, 7, 8, X, X, X, X}},
2049 {{5, 8, 9, 6, X, X, X, X}},
2050 {{X, X, X, X, X, X, X, X}},
2051 {{X, X, X, X, X, X, X, X}},
2052 }};
2053 return table_68;
2054 }
2055 case 1:
2056 {
2057 // new line is (5,7)
2058 static constexpr ndarray<unsigned int, 10, 8> table_57 = {{
2059 {{0, 4, 6, 7, X, X, X, X}},
2060 {{4, 1, 5, 8, X, X, X, X}},
2061 {{6, 5, 2, 9, X, X, X, X}},
2062 {{7, 8, 9, 3, X, X, X, X}},
2063 {{4, 5, 6, 7, X, X, X, X}},
2064 {{4, 7, 8, 5, X, X, X, X}},
2065 {{6, 9, 7, 5, X, X, X, X}},
2066 {{5, 8, 9, 7, X, X, X, X}},
2067 {{X, X, X, X, X, X, X, X}},
2068 {{X, X, X, X, X, X, X, X}},
2069 }};
2070 return table_57;
2071 }
2072 case 2:
2073 {
2074 // new line is (4,9)
2075 static constexpr ndarray<unsigned int, 10, 8> table_49 = {{
2076 {{0, 4, 6, 7, X, X, X, X}},
2077 {{4, 1, 5, 8, X, X, X, X}},
2078 {{6, 5, 2, 9, X, X, X, X}},
2079 {{7, 8, 9, 3, X, X, X, X}},
2080 {{4, 5, 6, 9, X, X, X, X}},
2081 {{4, 7, 8, 9, X, X, X, X}},
2082 {{6, 9, 7, 4, X, X, X, X}},
2083 {{5, 8, 9, 4, X, X, X, X}},
2084 {{X, X, X, X, X, X, X, X}},
2085 {{X, X, X, X, X, X, X, X}},
2086 }};
2087 return table_49;
2088 }
2089
2090 default:
2091 DEAL_II_ASSERT_UNREACHABLE();
2092 }
2093 case ReferenceCells::Pyramid:
2094 {
2095 static constexpr ndarray<unsigned int, 10, 8>
2096 cell_vertices_pyramid = {{
2097 {{0, 7, 5, 13, 9, X, X, X}}, // bottom pyramid
2098 {{7, 1, 13, 6, 10, X, X, X}}, //
2099 {{5, 13, 2, 8, 11, X, X, X}}, //
2100 {{13, 6, 8, 3, 12, X, X, X}}, //
2101 {{5, 13, 11, 9, X, X, X, X}}, // bottom wedges
2102 {{13, 6, 12, 10, X, X, X, X}}, //
2103 {{7, 13, 9, 10, X, X, X, X}}, //
2104 {{13, 8, 11, 12, X, X, X, X}}, //
2105 {{9, 11, 10, 12, 13, X, X, X}}, // upside down pyramid
2106 {{9, 10, 11, 12, 4, X, X, X}}, // top pyramid
2107 }};
2108 return cell_vertices_pyramid;
2109 }
2110
2111 case ReferenceCells::Wedge:
2112 {
2113 static constexpr ndarray<unsigned int, 10, 8>
2114 cell_vertices_wedge = {{
2115 {{0, 6, 8, 12, 15, 17, X, X}}, // bottom children
2116 {{6, 1, 7, 15, 13, 16, X, X}}, //
2117 {{8, 7, 2, 17, 16, 14, X, X}}, //
2118 {{7, 8, 6, 16, 17, 15, X, X}}, //
2119 {{12, 15, 17, 3, 9, 11, X, X}}, // top children
2120 {{15, 13, 16, 9, 4, 10, X, X}}, //
2121 {{17, 16, 14, 11, 10, 5, X, X}}, //
2122 {{16, 17, 15, 10, 11, 9, X, X}}, //
2123 {{X, X, X, X, X, X, X, X}},
2124 {{X, X, X, X, X, X, X, X}},
2125 }};
2126 return cell_vertices_wedge;
2127 }
2128
2129 case ReferenceCells::Hexahedron:
2130 {
2131 static constexpr ndarray<unsigned int, 10, 8> cell_vertices_hex =
2132 {{
2133 {{0, 10, 8, 24, 16, 22, 20, 26}}, // bottom children
2134 {{10, 1, 24, 9, 22, 17, 26, 21}}, //
2135 {{8, 24, 2, 11, 20, 26, 18, 23}}, //
2136 {{24, 9, 11, 3, 26, 21, 23, 19}}, //
2137 {{16, 22, 20, 26, 4, 14, 12, 25}}, // top children
2138 {{22, 17, 26, 21, 14, 5, 25, 13}}, //
2139 {{20, 26, 18, 23, 12, 25, 6, 15}}, //
2140 {{26, 21, 23, 19, 25, 13, 15, 7}}, //
2141 {{X, X, X, X, X, X, X, X}},
2142 {{X, X, X, X, X, X, X, X}},
2143 }};
2144 return cell_vertices_hex;
2145 }
2146
2147 default:
2148 DEAL_II_ASSERT_UNREACHABLE();
2149 }
2150 }
2151 else
2152 // We should never get here except in the 3d case.
2153 DEAL_II_ASSERT_UNREACHABLE();
2154
2155 static constexpr ndarray<unsigned int, 10, 8> empty{};
2156 return empty;
2157}
2158
2159
2160
2161template <int dim>
2162inline constexpr bool
2163ReferenceCell<dim>::is_hyper_cube() const
2164{
2165 return (*this == ReferenceCells::Vertex || *this == ReferenceCells::Line ||
2166 *this == ReferenceCells::Quadrilateral ||
2167 *this == ReferenceCells::Hexahedron);
2168}
2169
2170
2171
2172template <int dim>
2173inline constexpr bool
2174ReferenceCell<dim>::is_simplex() const
2175{
2176 return (*this == ReferenceCells::Vertex || *this == ReferenceCells::Line ||
2177 *this == ReferenceCells::Triangle ||
2178 *this == ReferenceCells::Tetrahedron);
2179}
2180
2181
2182
2183template <int dim>
2184inline constexpr unsigned int
2185ReferenceCell<dim>::get_dimension() const
2186{
2187 return dim;
2188}
2189
2190
2191
2192template <int dim>
2193Quadrature<dim>
2194ReferenceCell<dim>::get_midpoint_quadrature() const
2195{
2196 return Quadrature<dim>(std::vector<Point<dim>>({barycenter()}),
2197 std::vector<double>({volume()}));
2198}
2199
2200
2201
2202template <int dim>
2203inline constexpr unsigned int
2204ReferenceCell<dim>::n_vertices() const
2205{
2206 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2207
2208 if constexpr (dim == 0)
2209 return 1;
2210 else if constexpr (dim == 1)
2211 return 2;
2212 else if constexpr (dim == 2)
2213 {
2214 switch (this->kind)
2215 {
2216 case ReferenceCells::Triangle:
2217 return 3;
2218 case ReferenceCells::Quadrilateral:
2219 return 4;
2220 default:
2221 DEAL_II_NOT_IMPLEMENTED();
2222 }
2223 }
2224 else if constexpr (dim == 3)
2225 {
2226 switch (this->kind)
2227 {
2228 case ReferenceCells::Tetrahedron:
2229 return 4;
2230 case ReferenceCells::Pyramid:
2231 return 5;
2232 case ReferenceCells::Wedge:
2233 return 6;
2234 case ReferenceCells::Hexahedron:
2235 return 8;
2236 default:
2237 DEAL_II_NOT_IMPLEMENTED();
2238 }
2239 }
2240 else
2241 DEAL_II_NOT_IMPLEMENTED();
2242
2243 return numbers::invalid_unsigned_int;
2244}
2245
2246
2247
2248template <int dim>
2249inline constexpr unsigned int
2250ReferenceCell<dim>::n_lines() const
2251{
2252 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2253
2254 if constexpr (dim == 0)
2255 return 0;
2256 else if constexpr (dim == 1)
2257 return 1;
2258 else if constexpr (dim == 2)
2259 {
2260 switch (this->kind)
2261 {
2262 case ReferenceCells::Triangle:
2263 return 3;
2264 case ReferenceCells::Quadrilateral:
2265 return 4;
2266 default:
2267 DEAL_II_NOT_IMPLEMENTED();
2268 }
2269 }
2270 else if constexpr (dim == 3)
2271 {
2272 switch (this->kind)
2273 {
2274 case ReferenceCells::Tetrahedron:
2275 return 6;
2276 case ReferenceCells::Pyramid:
2277 return 8;
2278 case ReferenceCells::Wedge:
2279 return 9;
2280 case ReferenceCells::Hexahedron:
2281 return 12;
2282 default:
2283 DEAL_II_NOT_IMPLEMENTED();
2284 }
2285 }
2286 else
2287 DEAL_II_NOT_IMPLEMENTED();
2288
2289 return numbers::invalid_unsigned_int;
2290}
2291
2292
2293
2294template <int dim>
2295Point<dim>
2296ReferenceCell<dim>::vertex(const unsigned int v) const
2297{
2298 AssertIndexRange(v, n_vertices());
2299 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2300
2301 if constexpr (dim == 0)
2302 return Point<dim>();
2303 else if constexpr (dim == 1)
2304 {
2305 static const Point<dim> vertices[2] = {
2306 Point<dim>(), // the origin
2307 Point<dim>::unit_vector(0) // unit point along x-axis
2308 };
2309 return vertices[v];
2310 }
2311 else if constexpr (dim == 2)
2312 {
2313 switch (this->kind)
2314 {
2315 case ReferenceCells::Triangle:
2316 {
2317 static const Point<dim> vertices[3] = {
2318 Point<dim>(), // the origin
2319 Point<dim>::unit_vector(0), // unit point along x-axis
2320 Point<dim>::unit_vector(1) // unit point along y-axis
2321 };
2322 return vertices[v];
2323 }
2324 case ReferenceCells::Quadrilateral:
2325 {
2326 static const Point<dim> vertices[4] = {
2327 // First the two points on the x-axis
2328 Point<dim>(),
2329 Point<dim>::unit_vector(0),
2330 // Then these two points shifted in the y-direction
2331 Point<dim>() + Point<dim>::unit_vector(1),
2332 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(1)};
2333 return vertices[v];
2334 }
2335 }
2336 }
2337 else if constexpr (dim == 3)
2338 {
2339 switch (this->kind)
2340 {
2341 case ReferenceCells::Tetrahedron:
2342 {
2343 static const Point<dim> vertices[4] = {
2344 Point<dim>(), // the origin
2345 Point<dim>::unit_vector(0), // unit point along x-axis
2346 Point<dim>::unit_vector(1), // unit point along y-axis
2347 Point<dim>::unit_vector(2) // unit point along z-axis
2348 };
2349 return vertices[v];
2350 }
2351 case ReferenceCells::Pyramid:
2352 {
2353 static const Point<dim> vertices[5] = {
2354 Point<dim>{-1.0, -1.0, 0.0},
2355 Point<dim>{+1.0, -1.0, 0.0},
2356 Point<dim>{-1.0, +1.0, 0.0},
2357 Point<dim>{+1.0, +1.0, 0.0},
2358 Point<dim>{+0.0, +0.0, 1.0}};
2359 return vertices[v];
2360 }
2361 case ReferenceCells::Wedge:
2362 {
2363 static const Point<dim> vertices[6] = {
2364 // First the three points on the triangular base of the
2365 // wedge:
2366 Point<dim>(),
2367 Point<dim>::unit_vector(0),
2368 Point<dim>::unit_vector(1),
2369 // And now everything shifted in the z-direction again
2370 Point<dim>() + Point<dim>::unit_vector(2),
2371 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(2),
2372 Point<dim>::unit_vector(1) + Point<dim>::unit_vector(2)};
2373 return vertices[v];
2374 }
2375 case ReferenceCells::Hexahedron:
2376 {
2377 static const Point<dim> vertices[8] = {
2378 // First the two points on the x-axis
2379 Point<dim>(),
2380 Point<dim>::unit_vector(0),
2381 // Then these two points shifted in the y-direction
2382 Point<dim>() + Point<dim>::unit_vector(1),
2383 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(1),
2384 // And now all four points shifted in the z-direction
2385 Point<dim>() + Point<dim>::unit_vector(2),
2386 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(2),
2387 Point<dim>() + Point<dim>::unit_vector(1) +
2388 Point<dim>::unit_vector(2),
2389 Point<dim>::unit_vector(0) + Point<dim>::unit_vector(1) +
2390 Point<dim>::unit_vector(2)};
2391 return vertices[v];
2392 }
2393 }
2394 }
2395
2396 DEAL_II_NOT_IMPLEMENTED();
2397 return Point<dim>();
2398}
2399
2400
2401
2402template <int dim>
2403inline constexpr unsigned int
2404ReferenceCell<dim>::n_faces() const
2405{
2406 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2407
2408 if constexpr (dim == 0)
2409 return 0;
2410 else if constexpr (dim == 1)
2411 return 2;
2412 else if constexpr (dim == 2)
2413 {
2414 switch (this->kind)
2415 {
2416 case ReferenceCells::Triangle:
2417 return 3;
2418 case ReferenceCells::Quadrilateral:
2419 return 4;
2420 }
2421 }
2422 else if constexpr (dim == 3)
2423 {
2424 switch (this->kind)
2425 {
2426 case ReferenceCells::Tetrahedron:
2427 return 4;
2428 case ReferenceCells::Pyramid:
2429 return 5;
2430 case ReferenceCells::Wedge:
2431 return 5;
2432 case ReferenceCells::Hexahedron:
2433 return 6;
2434 }
2435 }
2436
2437 DEAL_II_NOT_IMPLEMENTED();
2438 return numbers::invalid_unsigned_int;
2439}
2440
2441
2442
2443template <int dim>
2444inline std_cxx20::ranges::iota_view<unsigned int, unsigned int>
2450
2451
2452
2453template <int dim>
2454inline std_cxx20::ranges::iota_view<unsigned int, unsigned int>
2455ReferenceCell<dim>::face_indices_by_type(
2456 const ReferenceCell<dim - 1> &face_ref_type) const
2457{
2458 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2459 Assert(face_ref_type != ReferenceCells::Invalid<dim - 1>,
2460 ExcNotImplemented());
2461
2462 if constexpr (dim == 0)
2463 {
2464 return this->face_indices(); // no faces
2465 }
2466 else if constexpr (dim == 1)
2467 {
2468 return this->face_indices();
2469 }
2470 else if constexpr (dim == 2)
2471 {
2472 return this->face_indices();
2473 }
2474 else if constexpr (dim == 3)
2475 {
2476 switch (this->kind)
2477 {
2478 case ReferenceCells::Tetrahedron:
2479 {
2480 if (face_ref_type == ReferenceCells::Triangle)
2481 return this->face_indices();
2482 else
2483 return std_cxx20::ranges::iota_view<unsigned int, unsigned int>(
2484 0U, 0U);
2485 }
2486 case ReferenceCells::Pyramid:
2487 {
2488 if (face_ref_type == ReferenceCells::Triangle)
2489 return std_cxx20::ranges::iota_view<unsigned int, unsigned int>(
2490 1U, 5U);
2491 else
2492 return std_cxx20::ranges::iota_view<unsigned int, unsigned int>(
2493 0U, 1U);
2494 }
2495 case ReferenceCells::Wedge:
2496 {
2497 if (face_ref_type == ReferenceCells::Triangle)
2498 return std_cxx20::ranges::iota_view<unsigned int, unsigned int>(
2499 0U, 2U);
2500 else
2501 return std_cxx20::ranges::iota_view<unsigned int, unsigned int>(
2502 2U, 5U);
2503 }
2504 case ReferenceCells::Hexahedron:
2505 {
2506 if (face_ref_type == ReferenceCells::Quadrilateral)
2507 return this->face_indices();
2508 else
2509 return std_cxx20::ranges::iota_view<unsigned int, unsigned int>(
2510 0U, 0U);
2511 }
2512 }
2513 }
2514
2515 DEAL_II_NOT_IMPLEMENTED();
2516 return std_cxx20::ranges::iota_view<unsigned int, unsigned int>(0U, 0U);
2517}
2518
2519
2520
2521template <int dim>
2522inline constexpr unsigned int
2523ReferenceCell<dim>::n_isotropic_refinement_choices() const
2524{
2525 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2526
2527 if constexpr (dim == 0)
2528 return 1;
2529 else if constexpr (dim == 1)
2530 return 1;
2531 else if constexpr (dim == 2)
2532 return 1;
2533 else if constexpr (dim == 3)
2534 {
2535 switch (this->kind)
2536 {
2537 case ReferenceCells::Tetrahedron:
2538 return 3;
2539 case ReferenceCells::Pyramid:
2540 return 1;
2541 case ReferenceCells::Wedge:
2542 return 1;
2543 case ReferenceCells::Hexahedron:
2544 return 1;
2545 }
2546 }
2547
2548 DEAL_II_NOT_IMPLEMENTED();
2549 return numbers::invalid_unsigned_int;
2550}
2551
2552
2553
2554template <int dim>
2555inline IsotropicRefinementChoice
2556ReferenceCell<dim>::get_isotropic_refinement_choice(
2557 const unsigned int ref_choice) const
2558{
2559 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2560 AssertIndexRange(ref_choice, n_isotropic_refinement_choices());
2561
2562 if constexpr (dim == 0)
2563 DEAL_II_NOT_IMPLEMENTED();
2564 else if constexpr (dim == 1)
2565 DEAL_II_NOT_IMPLEMENTED();
2566 else if constexpr (dim == 2)
2567 DEAL_II_NOT_IMPLEMENTED();
2568 else if constexpr (dim == 3)
2569 {
2570 switch (this->kind)
2571 {
2572 case ReferenceCells::Tetrahedron:
2573 {
2574 constexpr ndarray<IsotropicRefinementChoice, 3>
2575 isotropic_ref_choices = {
2576 {IsotropicRefinementChoice::cut_tet_68,
2577 IsotropicRefinementChoice::cut_tet_57,
2578 IsotropicRefinementChoice::cut_tet_49}};
2579 return isotropic_ref_choices[ref_choice];
2580 }
2581 default:
2582 DEAL_II_NOT_IMPLEMENTED();
2583 }
2584 }
2585
2586 DEAL_II_NOT_IMPLEMENTED();
2587 return IsotropicRefinementChoice::isotropic_refinement;
2588}
2589
2590
2591
2592template <int dim>
2593inline constexpr unsigned int
2594ReferenceCell<dim>::n_isotropic_children() const
2595{
2596 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2597
2598 if constexpr (dim == 0)
2599 return 1;
2600 else if constexpr (dim == 1)
2601 return 2;
2602 else if constexpr (dim == 2)
2603 return 4;
2604 else if constexpr (dim == 3)
2605 {
2606 switch (this->kind)
2607 {
2608 case ReferenceCells::Tetrahedron:
2609 return 8;
2610 case ReferenceCells::Pyramid:
2611 return 10;
2612 case ReferenceCells::Wedge:
2613 return 8;
2614 case ReferenceCells::Hexahedron:
2615 return 8;
2616 }
2617 }
2618
2619 DEAL_II_NOT_IMPLEMENTED();
2620 return numbers::invalid_unsigned_int;
2621}
2622
2623
2624
2625template <int dim>
2626unsigned int
2627ReferenceCell<dim>::n_children(const RefinementCase<dim> ref_case) const
2628{
2629 // Use GeometryInfo here to keep it the single source of truth
2630 if (this->is_hyper_cube())
2631 return GeometryInfo<dim>::n_children(ref_case);
2632
2633 // TODO: Rewrite once refinement cases are cleaned up
2634 // For all other cells, we assume that only isotropic refinement is
2635 // implemented. This is a strong assumption that should be removed as soon as
2636 // the issue regarding the uniqueness of refinement cases is resolved.
2637 //
2638 // This is a short outline of the issue: Different refinement options are
2639 // stored as an enum, primarily focusing on hypercubes. They therefore are (to
2640 // some extent) unsuitable for other cell types. For example, there are
2641 // multiple ways to isotropically refine a tetrahedron. Distinguishing between
2642 // these required the introduction of the additional
2643 // `IsotropicRefinementChoice` enum. However, this leads to a conflict where
2644 // the same integer value is assigned to multiple different refinement options
2645 // across the existing refinement cases. The most critical conflict is likely
2646 // `IsotropicRefinementChoice::isotropic_refinement` shadowing
2647 // `RefinementPossibilities<~>::no_refinement`.
2648 // The following code attempts to circumvent this mess by predicting whether
2649 // the given refinement case is isotropic and then returning either the number
2650 // of isotropic children or `0` (see below for why `0` is used).
2651
2652 // Detect isotropic refinement.
2653 // Only in 3D and only for Tetrahedrons we have ambiguity. But first assume
2654 // that we don't have a tet.
2655 bool is_isotropic =
2656 (ref_case == RefinementPossibilities<dim>::isotropic_refinement);
2657
2658 // check if assumption of line above was correct, else apply fix for tet
2659 if constexpr (dim == 3) // check for dim first (only 3D relevant for tet)
2660 if (*this == ReferenceCells::Tetrahedron) // tets
2661 is_isotropic =
2662 ((ref_case == //
2663 RefinementPossibilities<dim>::isotropic_refinement) ||
2664 // We must not check
2665 // `IsotropicRefinementChoice::isotropic_refinement`
2666 // because it coincides with
2667 // `RefinementPossibilities<dim>::no_refinement`
2668 // (both are 0).
2669 (ref_case ==
2670 static_cast<typename RefinementPossibilities<dim>::Possibilities>(
2671 IsotropicRefinementChoice::cut_tet_68)) || //
2672 (ref_case ==
2673 static_cast<typename RefinementPossibilities<dim>::Possibilities>(
2674 IsotropicRefinementChoice::cut_tet_57)) || //
2675 (ref_case ==
2676 static_cast<typename RefinementPossibilities<dim>::Possibilities>(
2677 IsotropicRefinementChoice::cut_tet_49)));
2678
2679
2680 // return the according value
2681 if (is_isotropic)
2682 return n_isotropic_children();
2683
2684 // Unsupported case
2685
2686 // DEAL_II_NOT_IMPLEMENTED();
2687 // We must return some data. Many parts of the code aren't aware
2688 // which refinement cases are actually supported, calling this
2689 // function with impossible refinement cases and expecting some
2690 // output. One example is the constructor of `FiniteElement`.
2691 return 0;
2692}
2693
2694
2695
2696template <int dim>
2697inline std_cxx20::ranges::iota_view<unsigned int, unsigned int>
2703
2704
2705
2706template <int dim>
2707inline std_cxx20::ranges::iota_view<unsigned int, unsigned int>
2713
2714
2715
2716template <int dim>
2717inline std_cxx20::ranges::iota_view<unsigned int, unsigned int>
2723
2724
2725
2726template <int dim>
2727inline ArrayView<const RefinementCase<dim>>
2728ReferenceCell<dim>::refinement_cases() const
2729{
2730 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2731
2732 if constexpr (dim == 0)
2733 {
2734 // Like equivalent_refinement_case(), to better enable generic
2735 // programming, return a value implying vertices can be refined since 1d
2736 // Triangulations can be refined.
2737 static constexpr std::array<RefinementCase<0>, 1> possibilities{
2738 {RefinementPossibilities<0>::isotropic_refinement}};
2739 return make_array_view(possibilities);
2740 }
2741 else if constexpr (dim == 1)
2742 {
2743 static constexpr std::array<RefinementCase<1>, 1> possibilities{
2744 {RefinementPossibilities<1>::isotropic_refinement}};
2745 return make_array_view(possibilities);
2746 }
2747 else if constexpr (dim == 2)
2748 {
2749 switch (this->kind)
2750 {
2751 case ReferenceCells::Triangle:
2752 {
2753 static constexpr std::array<RefinementCase<2>, 1> possibilities{
2754 {RefinementPossibilities<2>::isotropic_refinement}};
2755 return make_array_view(possibilities);
2756 }
2757 case ReferenceCells::Quadrilateral:
2758 {
2759 static constexpr std::array<RefinementCase<2>, 3> possibilities{
2760 {RefinementPossibilities<2>::cut_x,
2761 RefinementPossibilities<2>::cut_y,
2762 RefinementPossibilities<2>::isotropic_refinement}};
2763 return make_array_view(possibilities);
2764 }
2765 }
2766 }
2767 else if constexpr (dim == 3)
2768 {
2769 switch (this->kind)
2770 {
2771 case ReferenceCells::Tetrahedron:
2772 case ReferenceCells::Pyramid:
2773 case ReferenceCells::Wedge:
2774 {
2775 static constexpr std::array<RefinementCase<3>, 1> possibilities{
2776 {RefinementPossibilities<3>::isotropic_refinement}};
2777 return make_array_view(possibilities);
2778 }
2779 case ReferenceCells::Hexahedron:
2780 {
2781 static constexpr std::array<RefinementCase<3>, 7> possibilities{
2782 {RefinementPossibilities<3>::cut_x,
2783 RefinementPossibilities<3>::cut_y,
2784 RefinementPossibilities<3>::cut_xy,
2785 RefinementPossibilities<3>::cut_z,
2786 RefinementPossibilities<3>::cut_xz,
2787 RefinementPossibilities<3>::cut_yz,
2788 RefinementPossibilities<3>::isotropic_refinement}};
2789 return make_array_view(possibilities);
2790 }
2791 }
2792 }
2793
2794 DEAL_II_NOT_IMPLEMENTED();
2795 return ArrayView<RefinementCase<dim>>();
2796}
2797
2798
2799
2800template <int dim>
2801inline constexpr ReferenceCell<dim - 1>
2802ReferenceCell<dim>::face_reference_cell(const unsigned int face_no) const
2803{
2804 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2805 AssertIndexRange(face_no, n_faces());
2806
2807 if constexpr (dim == 0)
2808 return ReferenceCells::Invalid<dim - 1>;
2809 else if constexpr (dim == 1)
2810 return ReferenceCells::Vertex;
2811 else if constexpr (dim == 2)
2812 return ReferenceCells::Line;
2813 else if constexpr (dim == 3)
2814 {
2815 switch (this->kind)
2816 {
2817 case ReferenceCells::Tetrahedron:
2818 return ReferenceCells::Triangle;
2819 case ReferenceCells::Pyramid:
2820 if (face_no == 0)
2821 return ReferenceCells::Quadrilateral;
2822 else
2823 return ReferenceCells::Triangle;
2824 case ReferenceCells::Wedge:
2825 if (face_no > 1)
2826 return ReferenceCells::Quadrilateral;
2827 else
2828 return ReferenceCells::Triangle;
2829 case ReferenceCells::Hexahedron:
2830 return ReferenceCells::Quadrilateral;
2831 }
2832 }
2833
2834 DEAL_II_NOT_IMPLEMENTED();
2835 return ReferenceCells::Invalid<dim - 1>;
2836}
2837
2838
2839
2840template <int dim>
2841inline unsigned int
2842ReferenceCell<dim>::child_cell_on_face(
2843 const unsigned int face,
2844 const unsigned int subface,
2845 const types::geometric_orientation combined_orientation) const
2846{
2847 (void)combined_orientation;
2848
2849 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2850 AssertIndexRange(face, n_faces());
2851 AssertIndexRange(subface, face_reference_cell(face).n_isotropic_children());
2852
2853 if constexpr (dim == 0)
2854 DEAL_II_NOT_IMPLEMENTED();
2855 else if constexpr (dim == 1)
2856 DEAL_II_NOT_IMPLEMENTED();
2857 else if constexpr (dim == 2)
2858 {
2859 switch (this->kind)
2860 {
2861 case ReferenceCells::Triangle:
2862 {
2863 static constexpr ndarray<unsigned int, 3, 2> subcells = {
2864 {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
2865
2866 Assert(combined_orientation ==
2867 numbers::default_geometric_orientation ||
2868 combined_orientation ==
2869 numbers::reverse_line_orientation,
2870 ExcInternalError());
2871 return subcells[face][combined_orientation ==
2872 numbers::default_geometric_orientation ?
2873 subface :
2874 1 - subface];
2875 }
2876 case ReferenceCells::Quadrilateral:
2877 {
2878 const auto [face_orientation, face_rotation, face_flip] =
2879 internal::split_face_orientation(combined_orientation);
2880
2881 return GeometryInfo<2>::child_cell_on_face(
2882 RefinementCase<2>(
2883 RefinementPossibilities<2>::isotropic_refinement),
2884 face,
2885 subface,
2886 face_orientation,
2887 face_flip,
2888 face_rotation);
2889 }
2890 }
2891 }
2892 else if constexpr (dim == 3)
2893 {
2894 switch (this->kind)
2895 {
2896 case ReferenceCells::Tetrahedron:
2897 case ReferenceCells::Pyramid:
2898 case ReferenceCells::Wedge:
2899 {
2900 DEAL_II_NOT_IMPLEMENTED();
2901 break;
2902 }
2903 case ReferenceCells::Hexahedron:
2904 {
2905 const auto [face_orientation, face_rotation, face_flip] =
2906 internal::split_face_orientation(combined_orientation);
2907
2908 return GeometryInfo<3>::child_cell_on_face(
2909 RefinementCase<3>(
2910 RefinementPossibilities<3>::isotropic_refinement),
2911 face,
2912 subface,
2913 face_orientation,
2914 face_flip,
2915 face_rotation);
2916 }
2917 }
2918 }
2919
2920 DEAL_II_NOT_IMPLEMENTED();
2921 return numbers::invalid_unsigned_int;
2922}
2923
2924
2925
2926template <int dim>
2927inline std::array<unsigned int, 2>
2928ReferenceCell<dim>::standard_vertex_to_face_and_vertex_index(
2929 const unsigned int vertex) const
2930{
2931 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
2932 AssertIndexRange(vertex, n_vertices());
2933 // Work around a GCC warning at higher optimization levels by making all of
2934 // these tables the same size
2935 constexpr unsigned int X = numbers::invalid_unsigned_int;
2936
2937 if constexpr (dim == 0)
2938 DEAL_II_NOT_IMPLEMENTED();
2939 else if constexpr (dim == 1)
2940 {
2941 static constexpr ndarray<unsigned int, 2, 2> table = {
2942 {{{0, 0}}, {{1, 0}}}};
2943
2944 return table[vertex];
2945 }
2946 else if constexpr (dim == 2)
2947 {
2948 switch (this->kind)
2949 {
2950 case ReferenceCells::Triangle:
2951 {
2952 static constexpr ndarray<unsigned int, 6, 2> table = {
2953 {{{0, 0}}, {{0, 1}}, {{1, 1}}, {{X, X}}, {{X, X}}, {{X, X}}}};
2954
2955 return table[vertex];
2956 }
2957 case ReferenceCells::Quadrilateral:
2958 {
2959 return GeometryInfo<2>::standard_quad_vertex_to_line_vertex_index(
2960 vertex);
2961 }
2962 }
2963 }
2964 else if constexpr (dim == 3)
2965 {
2966 switch (this->kind)
2967 {
2968 case ReferenceCells::Tetrahedron:
2969 {
2970 static constexpr ndarray<unsigned int, 6, 2> table = {
2971 {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 2}}, {{X, X}}, {{X, X}}}};
2972
2973 return table[vertex];
2974 }
2975 case ReferenceCells::Pyramid:
2976 {
2977 static constexpr ndarray<unsigned int, 6, 2> table = {
2978 {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{0, 3}}, {{1, 2}}, {{X, X}}}};
2979
2980 return table[vertex];
2981 }
2982 case ReferenceCells::Wedge:
2983 {
2984 static constexpr ndarray<unsigned int, 6, 2> table = {
2985 {{{0, 1}}, {{0, 0}}, {{0, 2}}, {{1, 0}}, {{1, 1}}, {{1, 2}}}};
2986
2987 return table[vertex];
2988 }
2989 case ReferenceCells::Hexahedron:
2990 {
2991 return GeometryInfo<3>::standard_hex_vertex_to_quad_vertex_index(
2992 vertex);
2993 }
2994 }
2995 }
2996
2997 DEAL_II_NOT_IMPLEMENTED();
2998 return {};
2999}
3000
3001
3002
3003template <int dim>
3004inline std::array<unsigned int, 2>
3005ReferenceCell<dim>::standard_line_to_face_and_line_index(
3006 const unsigned int line) const
3007{
3008 AssertIndexRange(line, n_lines());
3009
3010 if constexpr (dim == 0)
3011 DEAL_II_NOT_IMPLEMENTED();
3012 else if constexpr (dim == 1)
3013 DEAL_II_NOT_IMPLEMENTED();
3014 else if constexpr (dim == 2)
3015 DEAL_II_NOT_IMPLEMENTED();
3016 else if constexpr (dim == 3)
3017 {
3018 switch (this->kind)
3019 {
3020 case ReferenceCells::Tetrahedron:
3021 {
3022 static const std::array<unsigned int, 2> table[6] = {
3023 {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}};
3024
3025 return table[line];
3026 }
3027 case ReferenceCells::Pyramid:
3028 {
3029 static const std::array<unsigned int, 2> table[8] = {{{0, 0}},
3030 {{0, 1}},
3031 {{0, 2}},
3032 {{0, 3}},
3033 {{1, 2}},
3034 {{2, 1}},
3035 {{1, 1}},
3036 {{2, 2}}};
3037
3038 return table[line];
3039 }
3040 case ReferenceCells::Wedge:
3041 {
3042 static const std::array<unsigned int, 2> table[9] = {{{0, 0}},
3043 {{0, 2}},
3044 {{0, 1}},
3045 {{1, 0}},
3046 {{1, 1}},
3047 {{1, 2}},
3048 {{2, 0}},
3049 {{2, 1}},
3050 {{3, 1}}};
3051
3052 return table[line];
3053 }
3054 case ReferenceCells::Hexahedron:
3055 {
3056 return GeometryInfo<3>::standard_hex_line_to_quad_line_index(
3057 line);
3058 }
3059 }
3060 }
3061
3062 DEAL_II_NOT_IMPLEMENTED();
3063 return {};
3064}
3065
3066
3067
3068template <int dim>
3069inline unsigned int
3070ReferenceCell<dim>::line_to_cell_vertices(const unsigned int line,
3071 const unsigned int vertex) const
3072{
3073 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3074 AssertIndexRange(vertex, 2);
3075 AssertIndexRange(line, n_lines());
3076
3077 if constexpr (dim == 0)
3078 return vertex;
3079 else if constexpr (dim == 1)
3080 return vertex;
3081 else if constexpr (dim == 2)
3082 {
3083 switch (this->kind)
3084 {
3085 case ReferenceCells::Triangle:
3086 {
3087 static constexpr ndarray<unsigned int, 3, 2> table = {
3088 {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
3089 return table[line][vertex];
3090 }
3091 case ReferenceCells::Quadrilateral:
3092 {
3093 static constexpr ndarray<unsigned int, 4, 2> table = {
3094 {{{0, 2}}, {{1, 3}}, {{0, 1}}, {{2, 3}}}};
3095 return table[line][vertex];
3096 }
3097 }
3098 }
3099 else if constexpr (dim == 3)
3100 {
3101 switch (this->kind)
3102 {
3103 case ReferenceCells::Tetrahedron:
3104 {
3105 static constexpr ndarray<unsigned int, 6, 2> table = {
3106 {{{0, 1}}, {{1, 2}}, {{2, 0}}, {{0, 3}}, {{1, 3}}, {{2, 3}}}};
3107 return table[line][vertex];
3108 }
3109 case ReferenceCells::Pyramid:
3110 {
3111 static constexpr ndarray<unsigned int, 8, 2> table = {{{{0, 2}},
3112 {{1, 3}},
3113 {{0, 1}},
3114 {{2, 3}},
3115 {{0, 4}},
3116 {{1, 4}},
3117 {{2, 4}},
3118 {{3, 4}}}};
3119 return table[line][vertex];
3120 }
3121 case ReferenceCells::Wedge:
3122 {
3123 static constexpr ndarray<unsigned int, 9, 2> table = {{{{0, 1}},
3124 {{1, 2}},
3125 {{2, 0}},
3126 {{3, 4}},
3127 {{4, 5}},
3128 {{5, 3}},
3129 {{0, 3}},
3130 {{1, 4}},
3131 {{2, 5}}}};
3132 return table[line][vertex];
3133 }
3134 case ReferenceCells::Hexahedron:
3135 {
3136 // first four lines comprise the bottom face, next four are the
3137 // top, and the last four are 'bottom to top'
3138 static constexpr ndarray<unsigned int, 12, 2> table = {
3139 {{{0, 2}},
3140 {{1, 3}},
3141 {{0, 1}},
3142 {{2, 3}},
3143 {{4, 6}},
3144 {{5, 7}},
3145 {{4, 5}},
3146 {{6, 7}},
3147 {{0, 4}},
3148 {{1, 5}},
3149 {{2, 6}},
3150 {{3, 7}}}};
3151 return table[line][vertex];
3152 }
3153 }
3154 }
3155
3156 DEAL_II_NOT_IMPLEMENTED();
3157 return numbers::invalid_unsigned_int;
3158}
3159
3160
3161
3162template <int dim>
3163inline unsigned int
3164ReferenceCell<dim>::face_to_cell_lines(
3165 const unsigned int face,
3166 const unsigned int line,
3167 const types::geometric_orientation combined_face_orientation) const
3168{
3169 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3170 AssertIndexRange(face, n_faces());
3171 AssertIndexRange(line, face_reference_cell(face).n_lines());
3172 (void)combined_face_orientation;
3173
3174 static constexpr unsigned int X = numbers::invalid_unsigned_int;
3175
3176 if constexpr (dim == 0)
3177 // We can't get here based on the assertions above: vertices
3178 // have no faces.
3179 DEAL_II_ASSERT_UNREACHABLE();
3180 else if constexpr (dim == 1)
3181 // We can't get here based on the assertions above: lines
3182 // have vertices as faces, which have no bounding lines.
3183 DEAL_II_ASSERT_UNREACHABLE();
3184 else if constexpr (dim == 2)
3185 // 2d cells have only one line per face, so the line we are
3186 // looking for is actually the face we were given:
3187 return face;
3188 else if constexpr (dim == 3)
3189 {
3190 switch (this->kind)
3191 {
3192 case ReferenceCells::Tetrahedron:
3193 {
3194 static constexpr ndarray<unsigned int, 4, 3> table = {
3195 {{{0, 1, 2}}, {{0, 3, 4}}, {{2, 5, 3}}, {{1, 4, 5}}}};
3196
3197 return table[face][standard_to_real_face_line(
3198 line, face, combined_face_orientation)];
3199 }
3200 case ReferenceCells::Pyramid:
3201 {
3202 static constexpr ndarray<unsigned int, 5, 4> table = {
3203 {{{0, 1, 2, 3}},
3204 {{0, 6, 4, X}},
3205 {{1, 5, 7, X}},
3206 {{2, 4, 5, X}},
3207 {{3, 7, 6, X}}}};
3208
3209 return table[face][standard_to_real_face_line(
3210 line, face, combined_face_orientation)];
3211 }
3212 case ReferenceCells::Wedge:
3213 {
3214 static constexpr ndarray<unsigned int, 5, 4> table = {
3215 {{{0, 2, 1, X}},
3216 {{3, 4, 5, X}},
3217 {{6, 7, 0, 3}},
3218 {{7, 8, 1, 4}},
3219 {{8, 6, 2, 5}}}};
3220
3221 return table[face][standard_to_real_face_line(
3222 line, face, combined_face_orientation)];
3223 }
3224 case ReferenceCells::Hexahedron:
3225 {
3226 const auto [face_orientation, face_rotation, face_flip] =
3227 internal::split_face_orientation(combined_face_orientation);
3228
3229 return GeometryInfo<3>::face_to_cell_lines(
3230 face, line, face_orientation, face_flip, face_rotation);
3231 }
3232 }
3233 }
3234
3235 DEAL_II_NOT_IMPLEMENTED();
3236 return numbers::invalid_unsigned_int;
3237}
3238
3239
3240
3241template <int dim>
3242inline unsigned int
3243ReferenceCell<dim>::face_to_cell_vertices(
3244 const unsigned int face,
3245 const unsigned int vertex,
3246 const types::geometric_orientation combined_face_orientation) const
3247{
3248 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3249 AssertIndexRange(face, n_faces());
3250 AssertIndexRange(vertex, face_reference_cell(face).n_vertices());
3251 AssertIndexRange(combined_face_orientation, n_face_orientations(face));
3252 AssertDimension(dim, get_dimension());
3253
3254 if constexpr (dim == 0)
3255 // We can't get here based on the assertions above: vertices
3256 // have no faces.
3257 DEAL_II_ASSERT_UNREACHABLE();
3258 else if constexpr (dim == 1)
3259 // A line has only one vertex per face, so the vertex we are
3260 // looking for is actually the face we were given:
3261 return face;
3262 else if constexpr (dim == 2)
3263 {
3264 switch (this->kind)
3265 {
3266 case ReferenceCells::Triangle:
3267 {
3268 static constexpr ndarray<unsigned int, 3, 2> table = {
3269 {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
3270
3271 return table[face][combined_face_orientation ==
3272 numbers::default_geometric_orientation ?
3273 vertex :
3274 (1 - vertex)];
3275 }
3276 case ReferenceCells::Quadrilateral:
3277 {
3278 const auto [face_orientation, face_rotation, face_flip] =
3279 internal::split_face_orientation(combined_face_orientation);
3280
3281 return GeometryInfo<2>::face_to_cell_vertices(
3282 face, vertex, face_orientation, face_flip, face_rotation);
3283 }
3284 }
3285 }
3286 else if constexpr (dim == 3)
3287 {
3288 switch (this->kind)
3289 {
3290 case ReferenceCells::Tetrahedron:
3291 {
3292 static constexpr ndarray<unsigned int, 4, 3> table = {
3293 {{{0, 1, 2}}, {{1, 0, 3}}, {{0, 2, 3}}, {{2, 1, 3}}}};
3294
3295 return table[face][standard_to_real_face_vertex(
3296 vertex, face, combined_face_orientation)];
3297 }
3298 case ReferenceCells::Pyramid:
3299 {
3300 constexpr auto X = numbers::invalid_unsigned_int;
3301 static constexpr ndarray<unsigned int, 5, 4> table = {
3302 {{{0, 1, 2, 3}},
3303 {{0, 2, 4, X}},
3304 {{3, 1, 4, X}},
3305 {{1, 0, 4, X}},
3306 {{2, 3, 4, X}}}};
3307
3308 return table[face][standard_to_real_face_vertex(
3309 vertex, face, combined_face_orientation)];
3310 }
3311 case ReferenceCells::Wedge:
3312 {
3313 constexpr auto X = numbers::invalid_unsigned_int;
3314 static constexpr ndarray<unsigned int, 6, 4> table = {
3315 {{{1, 0, 2, X}},
3316 {{3, 4, 5, X}},
3317 {{0, 1, 3, 4}},
3318 {{1, 2, 4, 5}},
3319 {{2, 0, 5, 3}}}};
3320
3321 return table[face][standard_to_real_face_vertex(
3322 vertex, face, combined_face_orientation)];
3323 }
3324 case ReferenceCells::Hexahedron:
3325 {
3326 const auto [face_orientation, face_rotation, face_flip] =
3327 internal::split_face_orientation(combined_face_orientation);
3328
3329 return GeometryInfo<3>::face_to_cell_vertices(
3330 face, vertex, face_orientation, face_flip, face_rotation);
3331 }
3332 }
3333 }
3334
3335 DEAL_II_NOT_IMPLEMENTED();
3336 return numbers::invalid_unsigned_int;
3337}
3338
3339
3340
3341template <int dim>
3342inline unsigned int
3343ReferenceCell<dim>::face_and_line_to_cell_vertices(
3344 const unsigned int face,
3345 const unsigned int line,
3346 const unsigned int vertex,
3347 const types::geometric_orientation combined_face_orientation) const
3348{
3349 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3350 AssertIndexRange(face, n_faces());
3351 AssertIndexRange(line, this->face_reference_cell(face).n_lines());
3352 AssertIndexRange(vertex, 2);
3353 AssertIndexRange(combined_face_orientation, n_face_orientations(face));
3354
3355 if constexpr (dim == 0)
3356 // We can't get here based on the assertions above: vertices
3357 // have no faces.
3358 DEAL_II_ASSERT_UNREACHABLE();
3359 else if constexpr (dim == 1)
3360 // We can't get here based on the assertions above: Lines'
3361 // faces are vertices that have no lines as sub-objects.
3362 DEAL_II_ASSERT_UNREACHABLE();
3363 else if constexpr (dim == 2)
3364 // For 2d objects, each face has (and is identical to) exactly
3365 // one line, so we can ask the other function:
3366 return face_to_cell_vertices(face, vertex, combined_face_orientation);
3367 else if constexpr (dim == 3)
3368 {
3369 switch (this->kind)
3370 {
3371 case ReferenceCells::Tetrahedron:
3372 {
3373 Assert(combined_face_orientation ==
3374 numbers::default_geometric_orientation,
3375 ExcNotImplemented());
3376
3377 const static ::ndarray<unsigned int, 4, 3, 2> table = {
3378 {{{{{0, 1}}, {{1, 2}}, {{2, 0}}}},
3379 {{{{1, 0}}, {{0, 3}}, {{3, 1}}}},
3380 {{{{0, 2}}, {{2, 3}}, {{3, 0}}}},
3381 {{{{2, 1}}, {{1, 3}}, {{3, 2}}}}}};
3382
3383 return table[face][line][vertex];
3384 }
3385 case ReferenceCells::Pyramid:
3386 {
3387 Assert(combined_face_orientation ==
3388 numbers::default_geometric_orientation,
3389 ExcNotImplemented());
3390
3391 static const unsigned int X = static_cast<unsigned int>(-1);
3392 const static ::ndarray<unsigned int, 5, 4, 2> table = {
3393 {{{{{0, 2}}, {{1, 3}}, {{0, 1}}, {{2, 3}}}},
3394 {{{{0, 2}}, {{2, 4}}, {{4, 0}}, {{X, X}}}},
3395 {{{{3, 1}}, {{1, 4}}, {{4, 3}}, {{X, X}}}},
3396 {{{{1, 0}}, {{0, 4}}, {{4, 1}}, {{X, X}}}},
3397 {{{{2, 3}}, {{3, 4}}, {{4, 2}}, {{X, X}}}}}};
3398
3399 return table[face][line][vertex];
3400 }
3401 case ReferenceCells::Wedge:
3402 {
3403 Assert(combined_face_orientation ==
3404 numbers::default_geometric_orientation,
3405 ExcNotImplemented());
3406
3407 static const unsigned int X = static_cast<unsigned int>(-1);
3408 const static ::ndarray<unsigned int, 5, 4, 2> table = {
3409 {{{{{1, 0}}, {{0, 2}}, {{2, 1}}, {{X, X}}}},
3410 {{{{3, 4}}, {{4, 5}}, {{5, 3}}, {{X, X}}}},
3411 {{{{0, 3}}, {{1, 4}}, {{0, 1}}, {{3, 4}}}},
3412 {{{{1, 4}}, {{2, 5}}, {{1, 2}}, {{4, 5}}}},
3413 {{{{2, 5}}, {{0, 3}}, {{2, 0}}, {{5, 3}}}}}};
3414
3415 return table[face][line][vertex];
3416 }
3417 case ReferenceCells::Hexahedron:
3418 {
3419 Assert(combined_face_orientation ==
3420 numbers::default_geometric_orientation,
3421 ExcNotImplemented());
3422
3423 const static ::ndarray<unsigned int, 6, 4, 2> table = {
3424 {{{{{0, 4}}, {{2, 6}}, {{0, 2}}, {{4, 6}}}},
3425 {{{{1, 5}}, {{3, 7}}, {{1, 3}}, {{5, 7}}}},
3426 {{{{0, 1}}, {{4, 5}}, {{0, 4}}, {{1, 5}}}},
3427 {{{{2, 3}}, {{6, 7}}, {{2, 6}}, {{3, 7}}}},
3428 {{{{0, 2}}, {{1, 3}}, {{0, 1}}, {{2, 3}}}},
3429 {{{{4, 6}}, {{5, 7}}, {{4, 5}}, {{6, 7}}}}}};
3430
3431 return table[face][line][vertex];
3432 }
3433 }
3434 }
3435
3436
3437 DEAL_II_NOT_IMPLEMENTED();
3438 return numbers::invalid_unsigned_int;
3439}
3440
3441
3442
3443template <int dim>
3444Point<dim>
3445ReferenceCell<dim>::face_vertex_location(const unsigned int face,
3446 const unsigned int vertex) const
3447{
3448 return this->vertex(face_to_cell_vertices(
3449 face, vertex, numbers::default_geometric_orientation));
3450}
3451
3452
3453
3454template <int dim>
3455inline unsigned int
3456ReferenceCell<dim>::standard_to_real_face_vertex(
3457 const unsigned int vertex,
3458 const unsigned int face,
3459 const types::geometric_orientation face_orientation) const
3460{
3461 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3462 AssertIndexRange(face, n_faces());
3463 AssertIndexRange(vertex, face_reference_cell(face).n_vertices());
3464 AssertIndexRange(face_orientation, n_face_orientations(face));
3465
3466 if constexpr (dim == 1)
3467 return vertex;
3468 else if constexpr (dim == 2)
3469 return line_vertex_permutations[face_orientation][vertex];
3470 else if constexpr (dim == 3)
3471 {
3472 switch (this->kind)
3473 {
3474 case ReferenceCells::Tetrahedron:
3475 return triangle_vertex_permutations[face_orientation][vertex];
3476
3477 case ReferenceCells::Pyramid:
3478 if (face == 0)
3479 return quadrilateral_vertex_permutations[face_orientation]
3480 [vertex];
3481 else
3482 return triangle_vertex_permutations[face_orientation][vertex];
3483
3484 case ReferenceCells::Wedge:
3485 if (face < 2)
3486 return triangle_vertex_permutations[face_orientation][vertex];
3487 else
3488 return quadrilateral_vertex_permutations[face_orientation]
3489 [vertex];
3490
3491 case ReferenceCells::Hexahedron:
3492 return quadrilateral_vertex_permutations[face_orientation][vertex];
3493 }
3494 }
3495
3496 DEAL_II_NOT_IMPLEMENTED();
3497 return numbers::invalid_unsigned_int;
3498}
3499
3500
3501
3502template <int dim>
3503inline unsigned int
3504ReferenceCell<dim>::standard_to_real_face_line(
3505 const unsigned int line,
3506 const unsigned int face,
3507 const types::geometric_orientation face_orientation) const
3508{
3509 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3510 AssertIndexRange(face, n_faces());
3511 AssertIndexRange(line, face_reference_cell(face).n_lines());
3512 AssertIndexRange(face_orientation, n_face_orientations(face));
3513
3514 if constexpr (dim == 0)
3515 DEAL_II_NOT_IMPLEMENTED();
3516 else if constexpr (dim == 1)
3517 DEAL_II_NOT_IMPLEMENTED();
3518 else if constexpr (dim == 2)
3519 DEAL_II_NOT_IMPLEMENTED();
3520 else if constexpr (dim == 3)
3521 {
3522 switch (this->kind)
3523 {
3524 case ReferenceCells::Tetrahedron:
3525 return triangle_line_permutations[face_orientation][line];
3526
3527 case ReferenceCells::Pyramid:
3528 if (face == 0)
3529 return quadrilateral_line_permutations[face_orientation][line];
3530 else
3531 return triangle_line_permutations[face_orientation][line];
3532
3533 case ReferenceCells::Wedge:
3534 if (face < 2)
3535 return triangle_line_permutations[face_orientation][line];
3536 else
3537 return quadrilateral_line_permutations[face_orientation][line];
3538
3539 case ReferenceCells::Hexahedron:
3540 return quadrilateral_line_permutations[face_orientation][line];
3541 }
3542 }
3543
3544 DEAL_II_NOT_IMPLEMENTED();
3545 return numbers::invalid_unsigned_int;
3546}
3547
3548
3549
3550namespace ReferenceCells
3551{
3552 template <int dim>
3553 constexpr const std::vector<ReferenceCell<dim>>
3554 get_reference_cells_in_dim()
3555 {
3556 if constexpr (dim == 0)
3557 return {{ReferenceCells::Vertex}};
3558 else if constexpr (dim == 1)
3559 return {{ReferenceCells::Line}};
3560 else if constexpr (dim == 2)
3561 return {{ReferenceCells::Triangle, ReferenceCells::Quadrilateral}};
3562 else if constexpr (dim == 3)
3563 return {{ReferenceCells::Tetrahedron,
3564 ReferenceCells::Pyramid,
3565 ReferenceCells::Wedge,
3566 ReferenceCells::Hexahedron}};
3567 else
3568 DEAL_II_NOT_IMPLEMENTED();
3569
3570 return {};
3571 }
3572
3573
3574
3575 template <int dim>
3576 inline constexpr const ReferenceCell<dim> &
3577 get_simplex()
3578 {
3579 if constexpr (dim == 0)
3580 return ReferenceCells::Vertex;
3581 else if constexpr (dim == 1)
3582 return ReferenceCells::Line;
3583 else if constexpr (dim == 2)
3584 return ReferenceCells::Triangle;
3585 else if constexpr (dim == 3)
3586 return ReferenceCells::Tetrahedron;
3587 else
3588 DEAL_II_NOT_IMPLEMENTED();
3589
3590 return ReferenceCells::Invalid<dim>;
3591 }
3592
3593
3594
3595 template <int dim>
3596 inline constexpr const ReferenceCell<dim> &
3597 get_hypercube()
3598 {
3599 if constexpr (dim == 0)
3600 return ReferenceCells::Vertex;
3601 else if constexpr (dim == 1)
3602 return ReferenceCells::Line;
3603 else if constexpr (dim == 2)
3604 return ReferenceCells::Quadrilateral;
3605 else if constexpr (dim == 3)
3606 return ReferenceCells::Hexahedron;
3607 else
3608 DEAL_II_NOT_IMPLEMENTED();
3609
3610 return ReferenceCells::Invalid<dim>;
3611 }
3612
3613
3614
3615 template <int dim>
3616 inline constexpr unsigned int
3617 n_reference_cells_in_dim()
3618 {
3619 if constexpr (dim == 0 || dim == 1)
3620 return 1; // Vertex or Line
3621 else if (dim == 2)
3622 return 2; // Quad and Tri
3623 else if (dim == 3)
3624 return 4; // Tet, Pyramid, Wedge and Hex
3625 else
3626 DEAL_II_NOT_IMPLEMENTED();
3627
3628 return 0;
3629 }
3630
3631
3632
3633 template <int dim>
3634 inline ReferenceCell<dim>
3635 n_vertices_to_reference_cell(const unsigned int n_vertices)
3636 {
3637 if constexpr (dim == 0)
3638 if (n_vertices == 1)
3639 return ReferenceCells::Vertex;
3640
3641 if constexpr (dim == 1)
3642 if (n_vertices == 2)
3643 return ReferenceCells::Line;
3644
3645 if constexpr (dim == 2)
3646 {
3647 if (n_vertices == 3)
3648 return ReferenceCells::Triangle;
3649 else if (n_vertices == 4)
3650 return ReferenceCells::Quadrilateral;
3651 }
3652
3653 if constexpr (dim == 3)
3654 {
3655 if (n_vertices == 4)
3656 return ReferenceCells::Tetrahedron;
3657 else if (n_vertices == 5)
3658 return ReferenceCells::Pyramid;
3659 if (n_vertices == 6)
3660 return ReferenceCells::Wedge;
3661 else if (n_vertices == 8)
3662 return ReferenceCells::Hexahedron;
3663 }
3664
3665 Assert(false,
3666 ExcMessage("The combination of dim = " + std::to_string(dim) +
3667 " and n_vertices = " + std::to_string(n_vertices) +
3668 " does not correspond to a known reference cell type."));
3669 return ReferenceCells::Invalid<dim>;
3670 }
3671
3672
3673
3674 template <int dim>
3675 inline constexpr unsigned int
3680
3681
3682
3683 template <int dim>
3684 inline constexpr unsigned int
3689
3690
3691
3692 template <int dim>
3693 inline constexpr unsigned int
3698
3699
3700
3701 template <int dim>
3702 inline constexpr unsigned int
3703 max_n_children()
3704 {
3705 // Pyramids (3D) have 10 children.
3706 return dim == 3 ? 10 : GeometryInfo<dim>::max_children_per_cell;
3707 }
3708} // namespace ReferenceCells
3709
3710
3711
3712template <int dim>
3713inline ReferenceCell<dim>
3714ReferenceCell<dim>::n_vertices_to_type(const int dim_,
3715 const unsigned int n_vertices)
3716{
3717 AssertDimension(dim, dim_);
3718 return ReferenceCells::n_vertices_to_reference_cell<dim>(n_vertices);
3719}
3720
3721
3722
3723template <int dim>
3724inline double
3725ReferenceCell<dim>::d_linear_shape_function(const Point<dim> &xi,
3726 const unsigned int i) const
3727{
3728 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3729 AssertIndexRange(i, n_vertices());
3730
3731 if constexpr (dim == 0)
3732 {
3733 (void)xi;
3734 return 1.0;
3735 }
3736 else if constexpr (dim == 1)
3737 {
3738 return GeometryInfo<dim>::d_linear_shape_function(xi, i);
3739 }
3740 else if constexpr (dim == 2)
3741 {
3742 switch (this->kind)
3743 {
3744 case ReferenceCells::Quadrilateral:
3745 return GeometryInfo<dim>::d_linear_shape_function(xi, i);
3746 case ReferenceCells::Triangle:
3747 {
3748 // see also BarycentricPolynomials<2>::compute_value
3749 switch (i)
3750 {
3751 case 0:
3752 return 1.0 - xi[std::min(0, dim - 1)] -
3753 xi[std::min(1, dim - 1)];
3754 case 1:
3755 return xi[std::min(0, dim - 1)];
3756 case 2:
3757 return xi[std::min(1, dim - 1)];
3758 default:
3759 DEAL_II_ASSERT_UNREACHABLE();
3760 }
3761 }
3762 }
3763 }
3764 else if constexpr (dim == 3)
3765 {
3766 switch (this->kind)
3767 {
3768 case ReferenceCells::Hexahedron:
3769 return GeometryInfo<dim>::d_linear_shape_function(xi, i);
3770
3771 case ReferenceCells::Tetrahedron:
3772 {
3773 // see also BarycentricPolynomials<3>::compute_value
3774 switch (i)
3775 {
3776 case 0:
3777 return 1.0 - xi[std::min(0, dim - 1)] -
3778 xi[std::min(1, dim - 1)] - xi[std::min(2, dim - 1)];
3779 case 1:
3780 return xi[std::min(0, dim - 1)];
3781 case 2:
3782 return xi[std::min(1, dim - 1)];
3783 case 3:
3784 return xi[std::min(2, dim - 1)];
3785 default:
3786 DEAL_II_ASSERT_UNREACHABLE();
3787 }
3788 }
3789
3790 case ReferenceCells::Pyramid:
3791 {
3792 // see also ScalarLagrangePolynomialPyramid::compute_value()
3793 const double Q14 = 0.25;
3794
3795 const double r = xi[std::min(0, dim - 1)];
3796 const double s = xi[std::min(1, dim - 1)];
3797 const double t = xi[std::min(2, dim - 1)];
3798
3799 const double ratio =
3800 (std::fabs(t - 1.0) > 1.0e-14 ? (r * s * t) / (1.0 - t) : 0.0);
3801
3802 if (i == 0)
3803 return Q14 * ((1.0 - r) * (1.0 - s) - t + ratio);
3804 if (i == 1)
3805 return Q14 * ((1.0 + r) * (1.0 - s) - t - ratio);
3806 if (i == 2)
3807 return Q14 * ((1.0 - r) * (1.0 + s) - t - ratio);
3808 if (i == 3)
3809 return Q14 * ((1.0 + r) * (1.0 + s) - t + ratio);
3810 else
3811 return t;
3812 }
3813
3814 case ReferenceCells::Wedge:
3815 // see also ScalarLagrangePolynomialWedge::compute_value()
3816 return ReferenceCells::Triangle.d_linear_shape_function(
3817 Point<2>(xi[std::min(0, dim - 1)],
3818 xi[std::min(1, dim - 1)]),
3819 i % 3) *
3820 ReferenceCells::Line.d_linear_shape_function(
3821 Point<1>(xi[std::min(2, dim - 1)]), i / 3);
3822 }
3823 }
3824
3825 DEAL_II_NOT_IMPLEMENTED();
3826 return 0.0;
3827}
3828
3829
3830
3831template <int dim>
3832inline Tensor<1, dim>
3833ReferenceCell<dim>::d_linear_shape_function_gradient(const Point<dim> &xi,
3834 const unsigned int i) const
3835{
3836 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3837
3838 if constexpr (dim == 0)
3839 return GeometryInfo<dim>::d_linear_shape_function_gradient(xi, i);
3840 else if constexpr (dim == 1)
3841 return GeometryInfo<dim>::d_linear_shape_function_gradient(xi, i);
3842 else if constexpr (dim == 2)
3843 {
3844 switch (this->kind)
3845 {
3846 case ReferenceCells::Quadrilateral:
3847 return GeometryInfo<dim>::d_linear_shape_function_gradient(xi, i);
3848
3849 case ReferenceCells::Triangle:
3850 // see also BarycentricPolynomials<2>::compute_grad()
3851 switch (i)
3852 {
3853 case 0:
3854 return Point<dim>(-1.0, -1.0);
3855 case 1:
3856 return Point<dim>(+1.0, +0.0);
3857 case 2:
3858 return Point<dim>(+0.0, +1.0);
3859 default:
3860 DEAL_II_ASSERT_UNREACHABLE();
3861 }
3862 }
3863 }
3864 else if constexpr (dim == 3)
3865 {
3866 switch (this->kind)
3867 {
3868 case ReferenceCells::Hexahedron:
3869 return GeometryInfo<dim>::d_linear_shape_function_gradient(xi, i);
3870 default:
3871 DEAL_II_NOT_IMPLEMENTED();
3872 }
3873 }
3874
3875 DEAL_II_NOT_IMPLEMENTED();
3876 return Point<dim>(+0.0, +0.0, +0.0);
3877}
3878
3879
3880
3881template <int dim>
3882inline double
3883ReferenceCell<dim>::volume() const
3884{
3885 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3886
3887 if constexpr (dim == 0)
3888 return 1;
3889 else if constexpr (dim == 1)
3890 return 1;
3891 else if constexpr (dim == 2)
3892 {
3893 switch (this->kind)
3894 {
3895 case ReferenceCells::Triangle:
3896 return 1. / 2.;
3897 case ReferenceCells::Quadrilateral:
3898 return 1;
3899 }
3900 }
3901 else if constexpr (dim == 3)
3902 {
3903 switch (this->kind)
3904 {
3905 case ReferenceCells::Tetrahedron:
3906 return 1. / 6.;
3907 case ReferenceCells::Pyramid:
3908 return 4. / 3.;
3909 case ReferenceCells::Wedge:
3910 return 1. / 2.;
3911 case ReferenceCells::Hexahedron:
3912 return 1;
3913 }
3914 }
3915
3916 DEAL_II_NOT_IMPLEMENTED();
3917 return 0.0;
3918}
3919
3920
3921
3922template <int dim>
3923inline double
3924ReferenceCell<dim>::face_measure(const unsigned int face_no) const
3925{
3926 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3927 AssertIndexRange(face_no, n_faces());
3928 constexpr double X = std::numeric_limits<double>::signaling_NaN();
3929
3930 if constexpr (dim == 0)
3931 return 1.0;
3932 else if constexpr (dim == 1)
3933 return 1.0;
3934 else if constexpr (dim == 2)
3935 {
3936 switch (this->kind)
3937 {
3938 case ReferenceCells::Quadrilateral:
3939 return 1.0;
3940 case ReferenceCells::Triangle:
3941 static const std::array<double, 5> tri_faces{
3942 {1.0, std::sqrt(2.0), 1.0, X, X}};
3943 return tri_faces[face_no];
3944 }
3945 }
3946 else if constexpr (dim == 3)
3947 {
3948 switch (this->kind)
3949 {
3950 case ReferenceCells::Hexahedron:
3951 return 1.0;
3952 case ReferenceCells::Tetrahedron:
3953 {
3954 static const std::array<double, 5> tet_faces{
3955 {0.5, 0.5, 0.5, std::sqrt(3.0) / 2.0, X}};
3956 return tet_faces[face_no];
3957 }
3958 case ReferenceCells::Pyramid:
3959 {
3960 static const std::array<double, 5> pyramid_faces{
3961 {4,
3962 std::sqrt(2.0),
3963 std::sqrt(2.0),
3964 std::sqrt(2.0),
3965 std::sqrt(2.0)}};
3966 return pyramid_faces[face_no];
3967 }
3968 case ReferenceCells::Wedge:
3969 {
3970 const std::array<double, 5> wedge_faces{
3971 {0.5, 0.5, 1.0, std::sqrt(2.0), 1.0}};
3972 return wedge_faces[face_no];
3973 }
3974 }
3975 }
3976
3977 DEAL_II_NOT_IMPLEMENTED();
3978 return X;
3979}
3980
3981
3982
3983template <int dim>
3984inline Point<dim>
3985ReferenceCell<dim>::barycenter() const
3986{
3987 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
3988
3989 if constexpr (dim == 0)
3990 return Point<dim>();
3991 else if constexpr (dim == 1)
3992 return Point<dim>(1. / 2.);
3993 else if constexpr (dim == 2)
3994 {
3995 switch (this->kind)
3996 {
3997 case ReferenceCells::Triangle:
3998 return Point<dim>(1. / 3., 1. / 3.);
3999 case ReferenceCells::Quadrilateral:
4000 return Point<dim>(1. / 2., 1. / 2.);
4001 }
4002 }
4003 else if constexpr (dim == 3)
4004 {
4005 switch (this->kind)
4006 {
4007 case ReferenceCells::Tetrahedron:
4008 return Point<dim>(1. / 4., 1. / 4., 1. / 4.);
4009 case ReferenceCells::Pyramid:
4010 return Point<dim>(0, 0, 1. / 4.);
4011 case ReferenceCells::Wedge:
4012 return Point<dim>(1. / 3, 1. / 3, 1. / 2.);
4013 case ReferenceCells::Hexahedron:
4014 return Point<dim>(1. / 2., 1. / 2., 1. / 2.);
4015 }
4016 }
4017
4018 DEAL_II_NOT_IMPLEMENTED();
4019 return Point<dim>();
4020}
4021
4022
4023
4024template <int dim>
4025inline bool
4026ReferenceCell<dim>::contains_point(const Point<dim> &p,
4027 const double tolerance) const
4028{
4029 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4030
4031 if constexpr (dim == 0)
4032 {
4033 // Vertices are special cases in that they do not actually
4034 // have coordinates. Error out if this function is called
4035 // with a vertex:
4036 Assert(false,
4037 ExcMessage("Vertices are zero-dimensional objects and "
4038 "as a consequence have no coordinates. You "
4039 "cannot meaningfully ask whether a point is "
4040 "inside a vertex (within a certain tolerance) "
4041 "without coordinate values."));
4042 return false;
4043 }
4044 else if constexpr (dim == 1)
4045 {
4046 for (unsigned int d = 0; d < dim; ++d)
4047 if ((p[d] < -tolerance) || (p[d] > 1 + tolerance))
4048 return false;
4049 return true;
4050 }
4051 else if constexpr (dim == 2)
4052 {
4053 switch (this->kind)
4054 {
4055 case ReferenceCells::Quadrilateral:
4056 {
4057 for (unsigned int d = 0; d < dim; ++d)
4058 if ((p[d] < -tolerance) || (p[d] > 1 + tolerance))
4059 return false;
4060 return true;
4061 }
4062 case ReferenceCells::Triangle:
4063 {
4064 // First make sure that we are in the first quadrant or octant
4065 for (unsigned int d = 0; d < dim; ++d)
4066 if (p[d] < -tolerance)
4067 return false;
4068
4069 // Now we also need to make sure that we are below the diagonal
4070 // line or plane that delineates the simplex. This diagonal is
4071 // given by sum(p[d])<=1, and a diagonal a distance eps away is
4072 // given by sum(p[d])<=1+eps*sqrt(d). (For example, the point at
4073 // (1,1) is a distance of 1/sqrt(2) away from the diagonal. That
4074 // is, its sum satisfies
4075 // sum(p[d]) = 2 <= 1 + (1/sqrt(2)) * sqrt(2)
4076 // in other words, it satisfies the predicate with eps=1/sqrt(2).)
4077 double sum = 0;
4078 for (unsigned int d = 0; d < dim; ++d)
4079 sum += p[d];
4080 return (sum <= 1 + tolerance * std::sqrt(1. * dim));
4081 }
4082 }
4083 }
4084 else if constexpr (dim == 3)
4085 {
4086 switch (this->kind)
4087 {
4088 case ReferenceCells::Hexahedron:
4089 {
4090 for (unsigned int d = 0; d < dim; ++d)
4091 if ((p[d] < -tolerance) || (p[d] > 1 + tolerance))
4092 return false;
4093 return true;
4094 }
4095
4096 case ReferenceCells::Tetrahedron:
4097 {
4098 // First make sure that we are in the first quadrant or octant
4099 for (unsigned int d = 0; d < dim; ++d)
4100 if (p[d] < -tolerance)
4101 return false;
4102
4103 // Now we also need to make sure that we are below the diagonal
4104 // line or plane that delineates the simplex. This diagonal is
4105 // given by sum(p[d])<=1, and a diagonal a distance eps away is
4106 // given by sum(p[d])<=1+eps*sqrt(d). (For example, the point at
4107 // (1,1) is a distance of 1/sqrt(2) away from the diagonal. That
4108 // is, its sum satisfies
4109 // sum(p[d]) = 2 <= 1 + (1/sqrt(2)) * sqrt(2)
4110 // in other words, it satisfies the predicate with eps=1/sqrt(2).)
4111 double sum = 0;
4112 for (unsigned int d = 0; d < dim; ++d)
4113 sum += p[d];
4114 return (sum <= 1 + tolerance * std::sqrt(1. * dim));
4115 }
4116
4117 case ReferenceCells::Pyramid:
4118 {
4119 // A pyramid only lives in the upper half-space:
4120 if (p[2] < -tolerance)
4121 return false;
4122
4123 // It also only lives in the space below z=1:
4124 if (p[2] > 1 + tolerance)
4125 return false;
4126
4127 // Within what's left of the space, a pyramid is a cone that
4128 // tapers towards the top. First compute the distance of the point
4129 // to the axis in the max norm (this is the right norm because the
4130 // vertices of the pyramid are at points +/-1, +/-1):
4131 const double distance_from_axis =
4132 std::max(std::fabs(p[0]), std::fabs(p[1]));
4133
4134 // We are inside the pyramid if the distance from the axis is less
4135 // than (1-z)
4136 return (distance_from_axis <= 1 + tolerance - p[2]);
4137 }
4138
4139 case ReferenceCells::Wedge:
4140 {
4141 // The wedge we use is a triangle extruded into the third
4142 // dimension by one unit. So we can use the same logic as for
4143 // triangles above (i.e., for the simplex above, using dim==2)
4144 // and then check the third dimension separately.
4145
4146 if ((p[0] < -tolerance) || (p[1] < -tolerance))
4147 return false;
4148
4149 const double sum = p[0] + p[1];
4150 if (sum > 1 + tolerance * std::sqrt(2.0))
4151 return false;
4152
4153 if (p[2] < -tolerance)
4154 return false;
4155 if (p[2] > 1 + tolerance)
4156 return false;
4157
4158 return true;
4159 }
4160 }
4161 }
4162
4163 DEAL_II_NOT_IMPLEMENTED();
4164 return false;
4165}
4166
4167
4168
4169template <int dim>
4170inline Tensor<1, dim>
4171ReferenceCell<dim>::face_tangent_vector(const unsigned int face_no,
4172 const unsigned int i) const
4173{
4174 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4175 AssertIndexRange(face_no, n_faces());
4176 AssertIndexRange(i, dim - 1);
4177
4178 if constexpr (dim == 0)
4179 // 0d cells (vertices) have no faces
4180 DEAL_II_ASSERT_UNREACHABLE();
4181 else if constexpr (dim == 1)
4182 // 1d cells (lines) have points as faces which have no
4183 // tangent space
4184 DEAL_II_ASSERT_UNREACHABLE();
4185 else if constexpr (dim == 2)
4186 {
4187 switch (this->kind)
4188 {
4189 case ReferenceCells::Quadrilateral:
4190 return GeometryInfo<dim>::unit_tangential_vectors[face_no][i];
4191 case ReferenceCells::Triangle:
4192 {
4193 constexpr std::array<Tensor<1, dim>, 3> table = {
4194 {Point<dim>(1, 0),
4195 Point<dim>(-numbers::SQRT1_2, +numbers::SQRT1_2),
4196 Point<dim>(0, -1)}};
4197
4198 return table[face_no];
4199 }
4200 }
4201 }
4202 else if constexpr (dim == 3)
4203 {
4204 switch (this->kind)
4205 {
4206 case ReferenceCells::Hexahedron:
4207 return GeometryInfo<dim>::unit_tangential_vectors[face_no][i];
4208 case ReferenceCells::Tetrahedron:
4209 {
4210 // We need std::pow(1.0/3.0, 0.25) in a constexpr context, but
4211 // that function isn't constexpr yet so hard-code the value and
4212 // assert that it is close
4213 constexpr double third_r4 = 0.7598356856515925473311877;
4214 Assert(std::abs((third_r4 * third_r4) * (third_r4 * third_r4) -
4215 1.0 / 3.0) < 1e-14,
4216 ExcInternalError());
4217 constexpr ndarray<Tensor<1, dim>, 4, 2> table = {
4218 {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
4219 {{Point<dim>(1, 0, 0), Point<dim>(0, 0, 1)}},
4220 {{Point<dim>(0, 0, 1), Point<dim>(0, 1, 0)}},
4221 {{Point<dim>(-third_r4, +third_r4, 0),
4222 Point<dim>(-third_r4, 0, +third_r4)}}}};
4223
4224 return table[face_no][i];
4225 }
4226 case ReferenceCells::Pyramid:
4227 {
4228 constexpr ndarray<Tensor<1, dim>, 5, 2> table = {
4229 {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
4230 {{Point<dim>(+numbers::SQRT1_2, 0, +numbers::SQRT1_2),
4231 Point<dim>(0, 1, 0)}},
4232 {{Point<dim>(+numbers::SQRT1_2, 0, -numbers::SQRT1_2),
4233 Point<dim>(0, 1, 0)}},
4234 {{Point<dim>(1, 0, 0),
4235 Point<dim>(0, +numbers::SQRT1_2, +numbers::SQRT1_2)}},
4236 {{Point<dim>(1, 0, 0),
4237 Point<dim>(0, +numbers::SQRT1_2, -numbers::SQRT1_2)}}}};
4238
4239 return table[face_no][i];
4240 }
4241 case ReferenceCells::Wedge:
4242 {
4243 constexpr ndarray<Tensor<1, dim>, 5, 2> table = {
4244 {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
4245 {{Point<dim>(1, 0, 0), Point<dim>(0, 1, 0)}},
4246 {{Point<dim>(1, 0, 0), Point<dim>(0, 0, 1)}},
4247 {{Point<dim>(-numbers::SQRT1_2, +numbers::SQRT1_2, 0),
4248 Point<dim>(0, 0, 1)}},
4249 {{Point<dim>(0, 0, 1), Point<dim>(0, 1, 0)}}}};
4250
4251 return table[face_no][i];
4252 }
4253 }
4254 }
4255
4256 DEAL_II_NOT_IMPLEMENTED();
4257 return {};
4258}
4259
4260
4261
4262template <int dim>
4263inline Tensor<1, dim>
4264ReferenceCell<dim>::face_normal_vector(const unsigned int face_no) const
4265{
4266 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4267 AssertIndexRange(face_no, n_faces());
4268
4269 if (is_hyper_cube())
4270 return GeometryInfo<dim>::unit_normal_vector[face_no];
4271 else if constexpr (dim == 2)
4272 {
4273 Assert(*this == ReferenceCells::Triangle, ExcInternalError());
4274
4275 // Return the rotated vector
4276 return cross_product_2d(face_tangent_vector(face_no, 0));
4277 }
4278 else if constexpr (dim == 3)
4279 return cross_product_3d(face_tangent_vector(face_no, 0),
4280 face_tangent_vector(face_no, 1));
4281
4282 DEAL_II_NOT_IMPLEMENTED();
4283 return {};
4284}
4285
4286
4287
4288template <int dim>
4289inline unsigned int
4290ReferenceCell<dim>::n_face_orientations(const unsigned int face_no) const
4291{
4292 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4293 AssertIndexRange(face_no, n_faces());
4294
4295 switch (face_reference_cell(face_no))
4296 {
4297 case ReferenceCells::Vertex:
4298 return 1;
4299 case ReferenceCells::Line:
4300 return 2;
4301 case ReferenceCells::Triangle:
4302 return 6;
4303 case ReferenceCells::Quadrilateral:
4304 return 8;
4305 default:
4306 DEAL_II_NOT_IMPLEMENTED();
4307 }
4308
4309 DEAL_II_ASSERT_UNREACHABLE();
4310 return numbers::invalid_unsigned_int;
4311}
4312
4313
4314
4315template <int dim>
4316inline types::geometric_orientation
4317ReferenceCell<dim>::face_to_cell_line_orientation(
4318 const unsigned int face_line_no,
4319 const unsigned int face_no,
4320 const types::geometric_orientation combined_face_orientation,
4321 const types::geometric_orientation line_orientation) const
4322{
4323 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4324 (void)face_line_no;
4325 (void)face_no;
4326 (void)combined_face_orientation;
4327 (void)line_orientation;
4328
4329 if constexpr (dim == 0)
4330 DEAL_II_NOT_IMPLEMENTED();
4331 else if constexpr (dim == 1)
4332 DEAL_II_NOT_IMPLEMENTED();
4333 else if constexpr (dim == 2)
4334 DEAL_II_NOT_IMPLEMENTED();
4335 else if constexpr (dim == 3)
4336 {
4337 constexpr auto D = numbers::default_geometric_orientation;
4338 constexpr auto R = numbers::reverse_line_orientation;
4339 constexpr unsigned int X = numbers::invalid_unsigned_int;
4340
4341 auto quad_face = [=]() {
4342 static constexpr ::ndarray<types::geometric_orientation, 2, 8>
4343 orientation_quad{
4344 {{{D, D, D, R, R, R, R, D}}, {{D, D, R, D, R, R, D, R}}}};
4345 // We use face_line_no / 2 here since lines i and i + 1 are parallel
4346 // and, on a given face, have the same relative orientations.
4347 const bool match =
4348 (line_orientation ==
4349 orientation_quad[face_line_no / 2][combined_face_orientation]);
4350
4351 return match ? numbers::default_geometric_orientation :
4352 numbers::reverse_line_orientation;
4353 };
4354
4355 constexpr ::ndarray<types::geometric_orientation, 2, 6>
4356 orientation_triangle{{{{D, R, D, R, D, R}}, {{R, D, R, D, R, D}}}};
4357
4358 switch (this->kind)
4359 {
4360 case ReferenceCells::Tetrahedron:
4361 {
4362 static constexpr ::ndarray<unsigned int, 4, 3>
4363 combined_lines{
4364 {{{0, 0, 0}}, {{X, 0, 1}}, {{X, 0, X}}, {{X, X, X}}}};
4365 const auto combined_line = combined_lines[face_no][face_line_no];
4366
4367 Assert(
4368 combined_line != X,
4369 ExcMessage(
4370 "This function can only be called for following face-line "
4371 "combinations: (0,0), (0,1), (0,2), (1,1), (1,2), (2,1)"));
4372
4373 const bool match =
4374 (line_orientation ==
4375 orientation_triangle[combined_line]
4376 [combined_face_orientation]);
4377 return match ? numbers::default_geometric_orientation :
4378 numbers::reverse_line_orientation;
4379 }
4380 case ReferenceCells::Wedge:
4381 if (face_no < 2)
4382 {
4383 // lines 0,1,2 (face 0, face_line_no 0,1,2) are not
4384 // oriented the same way in the cell and on the face
4385 // lines 4,5,6 (face 1, face_line 0,1,2) are
4386 // oriented the same way in the cell and on the faces
4387 static constexpr ::ndarray<unsigned int, 2, 3>
4388 combined_lines{{{{1, 1, 1}}, {{0, 0, 0}}}};
4389 const auto combined_line =
4390 combined_lines[face_no][face_line_no];
4391
4392 const bool match =
4393 (line_orientation ==
4394 orientation_triangle[combined_line]
4395 [combined_face_orientation]);
4396 return match ? numbers::default_geometric_orientation :
4397 numbers::reverse_line_orientation;
4398 }
4399 else
4400 return quad_face();
4401 case ReferenceCells::Pyramid:
4402 if (face_no == 0)
4403 return quad_face();
4404 else
4405 {
4406 // face_line_no 0 already done above
4407 Assert(face_line_no == 1 || face_line_no == 2,
4408 ExcMessage(
4409 "Invalid face line number for the triangular faces"));
4410
4411 // face 0 with lines 0-3 already done above
4412 // lines 6 (face 1, face_line_no 1) and 5 (face 2, face_line_no
4413 // 1) are oriented the same way in the cell and on the faces
4414 // lines 4 (face 1, face_line 2) and 7 (face 2, face_line 2) are
4415 // not oriented the same way in the cell and on the faces
4416 static constexpr ::ndarray<unsigned int, 3, 3>
4417 combined_lines{{{{X, X, X}}, {{X, 0, 1}}, {{X, 0, 1}}}};
4418
4419 const auto combined_line =
4420 combined_lines[face_no][face_line_no];
4421
4422 const bool match =
4423 (line_orientation ==
4424 orientation_triangle[combined_line]
4425 [combined_face_orientation]);
4426 return match ? numbers::default_geometric_orientation :
4427 numbers::reverse_line_orientation;
4428 }
4429 case ReferenceCells::Hexahedron:
4430 return quad_face();
4431 }
4432 }
4433
4434 DEAL_II_NOT_IMPLEMENTED();
4435 return numbers::invalid_geometric_orientation;
4436}
4437
4438
4439
4440template <int dim>
4441inline unsigned int
4442ReferenceCell<dim>::opposite_face_index(const unsigned int face_no) const
4443{
4444 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4445 AssertIndexRange(face_no, n_faces());
4446
4447 if constexpr (dim == 0)
4448 DEAL_II_NOT_IMPLEMENTED();
4449 else if constexpr (dim == 1)
4450 return (face_no + 1) % (this->get_dimension() + 1);
4451 else if constexpr (dim == 2)
4452 {
4453 switch (this->kind)
4454 {
4455 case ReferenceCells::Quadrilateral:
4456 return GeometryInfo<2>::opposite_face[face_no];
4457 case ReferenceCells::Triangle:
4458 return (face_no + 1) % (this->get_dimension() + 1);
4459 }
4460 }
4461 else if constexpr (dim == 3)
4462 {
4463 switch (this->kind)
4464 {
4465 case ReferenceCells::Hexahedron:
4466 return GeometryInfo<3>::opposite_face[face_no];
4467 case ReferenceCells::Tetrahedron:
4468 return (face_no + 1) % (this->get_dimension() + 1);
4469 case ReferenceCells::Pyramid:
4470 {
4471 switch (face_no)
4472 {
4473 case 0:
4474 return 0;
4475 case 1:
4476 return 2;
4477 case 2:
4478 return 1;
4479 case 3:
4480 return 4;
4481 case 4:
4482 return 3;
4483 default:
4484 DEAL_II_ASSERT_UNREACHABLE();
4485 return numbers::invalid_unsigned_int;
4486 }
4487 }
4488 case ReferenceCells::Wedge:
4489 {
4490 switch (face_no)
4491 {
4492 case 0:
4493 return 1;
4494 case 1:
4495 return 0;
4496 case 2:
4497 return 3;
4498 case 3:
4499 return 4;
4500 case 4:
4501 return 2;
4502 default:
4503 DEAL_II_ASSERT_UNREACHABLE();
4504 return numbers::invalid_unsigned_int;
4505 }
4506 }
4507 }
4508 }
4509
4510 DEAL_II_NOT_IMPLEMENTED();
4511 return numbers::invalid_unsigned_int;
4512}
4513
4514
4515
4516namespace internal
4517{
4518 template <int dim, typename T>
4519 class NoPermutation : public ::ExceptionBase
4520 {
4521 public:
4535
4539 virtual ~NoPermutation() noexcept override = default;
4540
4544 virtual void
4545 print_info(std::ostream &out) const override
4546 {
4547 out << '[';
4548
4549 const unsigned int n_vertices = entity_type.n_vertices();
4550
4551 for (unsigned int i = 0; i < n_vertices; ++i)
4552 {
4553 out << vertices_0[i];
4554 if (i + 1 != n_vertices)
4555 out << ',';
4556 }
4557
4558 out << "] is not a valid permutation of [";
4559
4560 for (unsigned int i = 0; i < n_vertices; ++i)
4561 {
4562 out << vertices_1[i];
4563 if (i + 1 != n_vertices)
4564 out << ',';
4565 }
4566
4567 out << "]." << std::endl;
4568 }
4569
4573 const ReferenceCell<dim> entity_type;
4574
4578 const ArrayView<const T> vertices_0;
4579
4583 const ArrayView<const T> vertices_1;
4584 };
4585
4586} // namespace internal
4587
4588
4589
4590template <int dim>
4591template <typename T, std::size_t N>
4592inline types::geometric_orientation
4593ReferenceCell<dim>::compute_orientation(
4594 const std::array<T, N> &vertices_0,
4595 const std::array<T, N> &vertices_1) const
4596{
4597 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4598 Assert(N >= n_vertices(),
4599 ExcMessage("The number of array elements must be equal to or "
4600 "greater than the number of vertices of the cell "
4601 "referenced by this object."));
4602
4603 // Call the non-deprecated function, taking care of calling it only with
4604 // those array elements that we actually care about (see the note
4605 // in the documentation about the arguments potentially being
4606 // larger arrays than necessary).
4607 return get_combined_orientation(
4608 make_array_view(vertices_0.begin(), vertices_0.begin() + n_vertices()),
4609 make_array_view(vertices_1.begin(), vertices_1.begin() + n_vertices()));
4610}
4611
4612
4613
4614template <int dim>
4615template <typename T>
4616types::geometric_orientation
4617ReferenceCell<dim>::get_combined_orientation(
4618 const ArrayView<const T> &vertices_0,
4619 const ArrayView<const T> &vertices_1) const
4620{
4621 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4622 Assert(vertices_0.size() == n_vertices(),
4623 ExcMessage("The number of array elements must be equal to "
4624 "the number of vertices of the cell "
4625 "referenced by this object."));
4626 Assert(vertices_1.size() == n_vertices(),
4627 ExcMessage("The number of array elements must be equal to "
4628 "the number of vertices of the cell "
4629 "referenced by this object."));
4630
4631 auto compute_orientation = [&](const auto &table,
4632 const unsigned int n_vertices) {
4633 for (types::geometric_orientation o = 0; o < table.size(); ++o)
4634 {
4635 bool match = true;
4636 for (unsigned int j = 0; j < n_vertices; ++j)
4637 {
4638 match &= (vertices_0[j] == vertices_1[table[o][j]]);
4639 if (match == false)
4640 break;
4641 }
4642
4643 if (match == true)
4644 return o;
4645 }
4646
4647 // Do not use `this` in Assert because nvcc when using C++20 assumes that
4648 // `this` is an integer and we get the following error: invalid type
4649 // argument of unary '*' (have 'int')
4650 [[maybe_unused]] const auto &ref_cell = *this;
4651 Assert(false,
4652 (internal::NoPermutation<dim, T>(ref_cell, vertices_0, vertices_1)));
4653 return std::numeric_limits<types::geometric_orientation>::max();
4654 };
4655
4656 if constexpr (dim == 0)
4657 return compute_orientation(vertex_vertex_permutations, 1);
4658 else if constexpr (dim == 1)
4659 return compute_orientation(line_vertex_permutations, 2);
4660 else if constexpr (dim == 2)
4661 {
4662 switch (this->kind)
4663 {
4664 case ReferenceCells::Triangle:
4665 return compute_orientation(triangle_vertex_permutations, 3);
4666 case ReferenceCells::Quadrilateral:
4667 return compute_orientation(quadrilateral_vertex_permutations, 4);
4668 }
4669 }
4670 else if constexpr (dim == 3)
4671 DEAL_II_NOT_IMPLEMENTED();
4672
4673 Assert(false,
4674 (internal::NoPermutation<dim, T>(*this, vertices_0, vertices_1)));
4675 return std::numeric_limits<types::geometric_orientation>::max();
4676}
4677
4678
4679
4680template <int dim>
4681template <typename T>
4682boost::container::small_vector<T, 8>
4683ReferenceCell<dim>::permute_by_combined_orientation(
4684 const ArrayView<const T> &vertices,
4685 const types::geometric_orientation orientation) const
4686{
4687 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4688 AssertDimension(vertices.size(), n_vertices());
4689
4690 const auto permute =
4691 [this, &vertices](const auto &table,
4692 const types::geometric_orientation orientation) {
4693 AssertIndexRange(orientation, table.size());
4694
4695 boost::container::small_vector<T, 8> result(this->n_vertices());
4696 for (unsigned int j = 0; j < result.size(); ++j)
4697 result[j] = vertices[table[orientation][j]];
4698 return result;
4699 };
4700
4701 if constexpr (dim == 0)
4702 return permute(vertex_vertex_permutations, orientation);
4703 else if constexpr (dim == 1)
4704 return permute(line_vertex_permutations, orientation);
4705 else if constexpr (dim == 2)
4706 {
4707 switch (this->kind)
4708 {
4709 case ReferenceCells::Triangle:
4710 return permute(triangle_vertex_permutations, orientation);
4711 case ReferenceCells::Quadrilateral:
4712 return permute(quadrilateral_vertex_permutations, orientation);
4713 }
4714 }
4715 else if constexpr (dim == 3)
4716 DEAL_II_NOT_IMPLEMENTED();
4717
4718 DEAL_II_NOT_IMPLEMENTED();
4719 return {};
4720}
4721
4722
4723
4724template <int dim>
4725inline types::geometric_orientation
4726ReferenceCell<dim>::get_inverse_combined_orientation(
4727 const types::geometric_orientation orientation) const
4728{
4729 Assert(*this != ReferenceCells::Invalid<dim>, ExcNotImplemented());
4730
4731 if constexpr (dim == 0)
4732 // Things are always default-oriented in 1D
4733 return orientation;
4734 else if constexpr (dim == 1)
4735 // the 1d orientations are the identity and a flip: i.e., the
4736 // identity and an involutory mapping
4737 return orientation;
4738 else if constexpr (dim == 2)
4739 {
4740 switch (this->kind)
4741 {
4742 case ReferenceCells::Triangle:
4743 {
4744 AssertIndexRange(orientation, 6);
4745 constexpr std::array<types::geometric_orientation, 6> inverses{
4746 {0, 1, 4, 3, 2, 5}};
4747 return inverses[orientation];
4748 }
4749 case ReferenceCells::Quadrilateral:
4750 {
4751 AssertIndexRange(orientation, 8);
4752 constexpr std::array<types::geometric_orientation, 8> inverses{
4753 {0, 1, 6, 3, 4, 5, 2, 7}};
4754 return inverses[orientation];
4755 }
4756 }
4757 }
4758 else
4759 {
4760 (void)orientation;
4761 DEAL_II_NOT_IMPLEMENTED();
4762 }
4763
4764 DEAL_II_NOT_IMPLEMENTED();
4765 return std::numeric_limits<types::geometric_orientation>::max();
4766}
4767
4768
4769DEAL_II_NAMESPACE_CLOSE
4770
4771#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:1002
ArrayView::size
std::size_t size() const
Definition array_view.h:737
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::unv_vertex_to_deal_vertex
unsigned int unv_vertex_to_deal_vertex(const unsigned int vertex_n) const
Definition reference_cell.cc:792
ReferenceCell::get_nodal_type_quadrature
const Quadrature< dim > & get_nodal_type_quadrature() const
Definition reference_cell.cc:564
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:3164
ReferenceCell::quadrilateral_vertex_permutations
static constexpr ndarray< unsigned int, 8, 4 > quadrilateral_vertex_permutations
Definition reference_cell.h:1180
ReferenceCell::is_simplex
constexpr bool is_simplex() const
Definition reference_cell.h:2174
ReferenceCell::volume
double volume() const
Definition reference_cell.h:3883
ReferenceCell::n_vertices
constexpr unsigned int n_vertices() const
Definition reference_cell.h:2204
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:4593
ReferenceCell::new_isotropic_child_cell_vertices
const ndarray< unsigned int, 10, 8 > & new_isotropic_child_cell_vertices(const unsigned int refinement_choice) const
Definition reference_cell.h:2025
ReferenceCell::get_isotropic_refinement_choice
IsotropicRefinementChoice get_isotropic_refinement_choice(const unsigned int ref_choice) const
Definition reference_cell.h:2556
ReferenceCell::face_to_cell_line_orientation
types::geometric_orientation face_to_cell_line_orientation(const unsigned int face_line_no, const unsigned int face_no, const types::geometric_orientation face_orientation, const types::geometric_orientation line_orientation) const
Convert a line orientation defined relative to a face to the canonical per-cell line orientation.
Definition reference_cell.h:4317
ReferenceCell::serialize
void serialize(Archive &archive, const unsigned int)
Definition reference_cell.h:1483
ReferenceCell::face_measure
double face_measure(const unsigned int face_no) const
Definition reference_cell.h:3924
ReferenceCell::is_hyper_cube
constexpr bool is_hyper_cube() const
Definition reference_cell.h:2163
ReferenceCell::line_indices
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
Definition reference_cell.h:2718
ReferenceCell::faces_for_given_vertex
ArrayView< const unsigned int > faces_for_given_vertex(const unsigned int vertex_index) const
Definition reference_cell.h:1501
ReferenceCell::operator!=
constexpr bool operator!=(const ReferenceCell< dim2 > &type) const
Definition reference_cell.h:1307
ReferenceCell::ucd_vertex_to_deal_vertex
unsigned int ucd_vertex_to_deal_vertex(const unsigned int vertex_n) const
Definition reference_cell.cc:732
ReferenceCell::vtk_quadratic_type
unsigned int vtk_quadratic_type() const
Definition reference_cell.cc:872
ReferenceCell::closest_point
Point< dim > closest_point(const Point< dim > &p) const
Definition reference_cell.cc:1439
ReferenceCell::get_midpoint_quadrature
Quadrature< dim > get_midpoint_quadrature() const
Definition reference_cell.h:2194
ReferenceCell::d_linear_shape_function
double d_linear_shape_function(const Point< dim > &xi, const unsigned int i) const
Definition reference_cell.h:3725
ReferenceCell::get_dimension
constexpr unsigned int get_dimension() const
Definition reference_cell.h:2185
ReferenceCell::n_vertices_to_type
static ReferenceCell n_vertices_to_type(const int dim_, const unsigned int n_vertices)
Definition reference_cell.h:3714
ReferenceCell::vtk_lagrange_type
unsigned int vtk_lagrange_type() const
Definition reference_cell.cc:911
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:4617
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:2928
ReferenceCell::refinement_cases
ArrayView< const RefinementCase< dim > > refinement_cases() const
Definition reference_cell.h:2728
ReferenceCell::get_gauss_type_quadrature
Quadrature< dim > get_gauss_type_quadrature(const unsigned n_points_1d) const
Definition reference_cell.cc:541
ReferenceCell::line_to_cell_vertices
unsigned int line_to_cell_vertices(const unsigned int line, const unsigned int vertex) const
Definition reference_cell.h:3070
ReferenceCell::kind
std::uint8_t kind
Definition reference_cell.h:1139
ReferenceCell::vertex_indices
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
Definition reference_cell.h:2708
ReferenceCell::face_reference_cell
constexpr ReferenceCell< dim - 1 > face_reference_cell(const unsigned int face_index) const
Definition reference_cell.h:2802
ReferenceCell::operator>>
friend std::istream & operator>>(std::istream &in, ReferenceCell< dim_ > &reference_cell)
ReferenceCell::n_children
unsigned int n_children(const RefinementCase< dim > ref_case=RefinementCase< dim >::isotropic_refinement) const
Definition reference_cell.h:2627
ReferenceCell::n_faces
constexpr unsigned int n_faces() const
Definition reference_cell.h:2404
ReferenceCell::ReferenceCell
constexpr ReferenceCell()
Definition reference_cell.h:1474
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:3504
ReferenceCell::barycenter
Point< dim > barycenter() const
Definition reference_cell.h:3985
ReferenceCell::n_isotropic_refinement_choices
constexpr unsigned int n_isotropic_refinement_choices() const
Definition reference_cell.h:2523
ReferenceCell::quadrilateral_line_permutations
static constexpr ndarray< unsigned int, 8, 4 > quadrilateral_line_permutations
Definition reference_cell.h:1206
ReferenceCell::vtk_vertex_to_deal_vertex
unsigned int vtk_vertex_to_deal_vertex(const unsigned int vertex_index) const
Definition reference_cell.cc:1145
ReferenceCell::line_vertex_permutations
static constexpr ndarray< unsigned int, 3, 4 > line_vertex_permutations
Definition reference_cell.h:1160
ReferenceCell::n_lines
constexpr unsigned int n_lines() const
Definition reference_cell.h:2250
ReferenceCell::vertex
Point< dim > vertex(const unsigned int v) const
Definition reference_cell.h:2296
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:3456
ReferenceCell::equivalent_refinement_case
std::pair< unsigned int, RefinementCase< dim - 1 > > equivalent_refinement_case(const types::geometric_orientation combined_face_orientation, const internal::SubfaceCase< dim > subface_case, const unsigned int subface_no) const
Definition reference_cell.cc:303
ReferenceCell::triangle_line_permutations
static constexpr ndarray< unsigned int, 6, 3 > triangle_line_permutations
Definition reference_cell.h:1194
ReferenceCell::opposite_face_index
unsigned int opposite_face_index(const unsigned int face_no) const
Return the index of the face opposite to a given face on the reference cell.
Definition reference_cell.h:4442
ReferenceCell::get_default_mapping
std::unique_ptr< Mapping< dim, spacedim > > get_default_mapping(const unsigned int degree) const
Definition reference_cell.cc:477
ReferenceCell::vertex_vertex_permutations
static constexpr ndarray< unsigned int, 1, 1 > vertex_vertex_permutations
Definition reference_cell.h:1154
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:3833
ReferenceCell::to_string
std::string to_string() const
Definition reference_cell.cc:82
ReferenceCell::n_face_orientations
unsigned int n_face_orientations(const unsigned int face_no) const
Definition reference_cell.h:4290
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
Definition reference_cell.cc:950
ReferenceCell::operator==
constexpr bool operator==(const ReferenceCell< dim2 > &type) const
Definition reference_cell.h:1288
ReferenceCell::triangle_vertex_permutations
static constexpr ndarray< unsigned int, 6, 4 > triangle_vertex_permutations
Definition reference_cell.h:1167
ReferenceCell::gmsh_element_type
unsigned int gmsh_element_type() const
Definition reference_cell.cc:1205
ReferenceCell::operator!=
constexpr bool operator!=(const ReferenceCell< dim > &type) const
Definition reference_cell.h:1297
ReferenceCell::isotropic_child_indices
std_cxx20::ranges::iota_view< unsigned int, unsigned int > isotropic_child_indices() const
Definition reference_cell.h:2698
ReferenceCell::operator<<
friend std::ostream & operator<<(std::ostream &out, const ReferenceCell< dim_ > &reference_cell)
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:3005
ReferenceCell::face_vertex_location
Point< dim > face_vertex_location(const unsigned int face, const unsigned int vertex) const
Definition reference_cell.h:3445
ReferenceCell::face_tangent_vector
Tensor< 1, dim > face_tangent_vector(const unsigned int face_no, const unsigned int i) const
Definition reference_cell.h:4171
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:3243
ReferenceCell::exodusii_vertex_to_deal_vertex
unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const
Definition reference_cell.cc:612
ReferenceCell::face_indices_by_type
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices_by_type(const ReferenceCell< dim - 1 > &face_ref_type) const
Definition reference_cell.h:2455
ReferenceCell::child_cell_on_face
unsigned int child_cell_on_face(const unsigned int face, const unsigned int subface, const types::geometric_orientation combined_orientation) const
Definition reference_cell.h:2842
ReferenceCell::get_default_linear_mapping
const Mapping< dim, spacedim > & get_default_linear_mapping() const
Definition reference_cell.cc:507
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:4683
ReferenceCell::get_inverse_combined_orientation
types::geometric_orientation get_inverse_combined_orientation(const types::geometric_orientation orientation) const
Definition reference_cell.h:4726
ReferenceCell::memory_consumption
static constexpr std::size_t memory_consumption()
Definition reference_cell.h:1492
ReferenceCell::face_normal_vector
Tensor< 1, dim > face_normal_vector(const unsigned int face_no) const
Definition reference_cell.h:4264
ReferenceCell::n_isotropic_children
constexpr unsigned int n_isotropic_children() const
Definition reference_cell.h:2594
ReferenceCell::subface_vertex_location
Point< dim > subface_vertex_location(const unsigned int face_no, const unsigned int subface_no, const unsigned int subface_vertex_no, const RefinementCase< dim - 1 > face_refinement_case) const
Definition reference_cell.cc:267
ReferenceCell::contains_point
bool contains_point(const Point< dim > &p, const double tolerance=0) const
Definition reference_cell.h:4026
ReferenceCell::exodusii_face_to_deal_face
unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const
Definition reference_cell.cc:671
ReferenceCell::face_and_line_to_cell_vertices
unsigned int face_and_line_to_cell_vertices(const unsigned int face, const unsigned int line, const unsigned int vertex, const types::geometric_orientation face_orientation) const
Definition reference_cell.h:3343
ReferenceCell::vtk_linear_type
unsigned int vtk_linear_type() const
Definition reference_cell.cc:833
ReferenceCell::operator==
constexpr bool operator==(const ReferenceCell< dim > &type) const
Definition reference_cell.h:1278
ReferenceCell::ReferenceCell
constexpr ReferenceCell(const std::uint8_t kind)
Definition reference_cell.h:1262
ReferenceCell::face_indices
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
Definition reference_cell.h:2445
ReferenceCell::new_isotropic_child_face_lines
const ndarray< unsigned int, 13, 4 > & new_isotropic_child_face_lines(const unsigned int refinement_choice) const
Definition reference_cell.h:1576
ReferenceCell::new_isotropic_child_cell_faces
const ndarray< unsigned int, 10, 6 > & new_isotropic_child_cell_faces(const unsigned int refinement_choice) const
Definition reference_cell.h:1882
ReferenceCell::new_isotropic_child_face_line_vertices
const ndarray< unsigned int, 13, 4, 2 > & new_isotropic_child_face_line_vertices(const unsigned int refinement_choice) const
Definition reference_cell.h:1724
internal::NoPermutation::vertices_1
const ArrayView< const T > vertices_1
Definition reference_cell.h:4583
internal::NoPermutation::print_info
virtual void print_info(std::ostream &out) const override
Definition reference_cell.h:4545
internal::NoPermutation::entity_type
const ReferenceCell< dim > entity_type
Definition reference_cell.h:4573
internal::NoPermutation::NoPermutation
NoPermutation(const ReferenceCell< dim > &entity_type, const ArrayView< const T > &vertices_0, const ArrayView< const T > &vertices_1)
Definition reference_cell.h:4525
internal::NoPermutation::~NoPermutation
virtual ~NoPermutation() noexcept override=default
internal::NoPermutation::vertices_0
const ArrayView< const T > vertices_0
Definition reference_cell.h:4578
DEAL_II_DEPRECATED
#define DEAL_II_DEPRECATED
Definition config.h:294
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:38
DEAL_II_DEPRECATED_WITH_COMMENT
#define DEAL_II_DEPRECATED_WITH_COMMENT(comment)
Definition config.h:295
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:39
DEAL_II_ASSERT_UNREACHABLE
#define DEAL_II_ASSERT_UNREACHABLE()
Definition exception_macros.h:781
DEAL_II_NOT_IMPLEMENTED
#define DEAL_II_NOT_IMPLEMENTED()
Definition exception_macros.h:705
geometry_info.h
StandardExceptions::ExcNotImplemented
static ::ExceptionBase & ExcNotImplemented()
Assert
#define Assert(cond, exc)
Definition exception_macros.h:543
AssertDimension
#define AssertDimension(dim1, dim2)
Definition exception_macros.h:829
AssertIndexRange
#define AssertIndexRange(index, range)
Definition exception_macros.h:901
StandardExceptions::ExcInternalError
static ::ExceptionBase & ExcInternalError()
StandardExceptions::ExcMessage
static ::ExceptionBase & ExcMessage(std::string arg1)
ReferenceCells::max_n_children
constexpr unsigned int max_n_children()
Definition reference_cell.h:3703
ReferenceCells::Hexahedron
constexpr ReferenceCell< 3 > Hexahedron
Definition reference_cell.h:1358
ReferenceCells::Quadrilateral
constexpr ReferenceCell< 2 > Quadrilateral
Definition reference_cell.h:1350
ReferenceCells::n_vertices_to_reference_cell
ReferenceCell< dim > n_vertices_to_reference_cell(const unsigned int n_vertices)
Definition reference_cell.h:3635
ReferenceCells::Line
constexpr ReferenceCell< 1 > Line
Definition reference_cell.h:1346
ReferenceCells::max_n_lines
constexpr unsigned int max_n_lines()
Definition reference_cell.h:3685
ReferenceCells::get_reference_cells_in_dim
constexpr const std::vector< ReferenceCell< dim > > get_reference_cells_in_dim()
Definition reference_cell.h:3554
ReferenceCells::Triangle
constexpr ReferenceCell< 2 > Triangle
Definition reference_cell.h:1348
ReferenceCells::max_n_vertices
constexpr unsigned int max_n_vertices()
Definition reference_cell.h:3676
ReferenceCells::Tetrahedron
constexpr ReferenceCell< 3 > Tetrahedron
Definition reference_cell.h:1352
ReferenceCells::get_simplex
constexpr const ReferenceCell< dim > & get_simplex()
Definition reference_cell.h:3577
ReferenceCells::max_n_faces
constexpr unsigned int max_n_faces()
Definition reference_cell.h:3694
ReferenceCells::Invalid
constexpr ReferenceCell< dim > Invalid
Definition reference_cell.h:1362
ReferenceCells::Pyramid
constexpr ReferenceCell< 3 > Pyramid
Definition reference_cell.h:1354
ReferenceCells::get_hypercube
constexpr const ReferenceCell< dim > & get_hypercube()
Definition reference_cell.h:3597
ReferenceCells::n_reference_cells_in_dim
constexpr unsigned int n_reference_cells_in_dim()
Definition reference_cell.h:3617
ReferenceCells::Wedge
constexpr ReferenceCell< 3 > Wedge
Definition reference_cell.h:1356
ReferenceCells::Vertex
constexpr ReferenceCell< 0 > Vertex
Definition reference_cell.h:1344
internal::split_face_orientation
std::tuple< bool, bool, bool > split_face_orientation(const types::geometric_orientation combined_orientation)
Definition tria_orientation.h:44
internal::make_reference_cell_from_int
constexpr ReferenceCell< dim > make_reference_cell_from_int(const std::uint8_t kind)
Definition reference_cell.h:1318
numbers::invalid_geometric_orientation
constexpr types::geometric_orientation invalid_geometric_orientation
Definition types.h:369
numbers::invalid_unsigned_int
constexpr unsigned int invalid_unsigned_int
Definition types.h:228
numbers::SQRT1_2
constexpr double SQRT1_2
Definition numbers.h:260
numbers::reverse_line_orientation
constexpr types::geometric_orientation reverse_line_orientation
Definition types.h:355
numbers::default_geometric_orientation
constexpr types::geometric_orientation default_geometric_orientation
Definition types.h:342
std_cxx20::ranges::iota_view
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:43
std
STL namespace.
std::min
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
Definition vectorization.h:6596
std::max
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
Definition vectorization.h:6579
std::sqrt
::VectorizedArray< Number, width > sqrt(const ::VectorizedArray< Number, width > &)
Definition vectorization.h:6505
std::abs
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)
Definition vectorization.h:6563
types::geometric_orientation
std::uint8_t geometric_orientation
Definition types.h:38
ndarray
typename internal::ndarray::HelperArray< T, Ns... >::type ndarray
Definition ndarray.h:105
IsotropicRefinementChoice
IsotropicRefinementChoice
Definition reference_cell.h:76
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< dim > &reference_cell)
Definition reference_cell.cc:1659
operator<<
std::ostream & operator<<(std::ostream &out, const ReferenceCell< dim > &reference_cell)
Definition reference_cell.cc:1645
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_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