Reference documentation for deal.II version GIT 5ac9d67d2d 2023-06-07 21:45: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\}}\)
reference_cell.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2020 - 2022 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_tria_reference_cell_h
17 #define dealii_tria_reference_cell_h
18 
19 #include <deal.II/base/config.h>
20 
23 #include <deal.II/base/ndarray.h>
24 #include <deal.II/base/point.h>
25 #include <deal.II/base/tensor.h>
26 #include <deal.II/base/utilities.h>
27 
28 #include <boost/container/small_vector.hpp>
29 
30 #include <iosfwd>
31 #include <string>
32 
34 
35 // Forward declarations
36 #ifndef DOXYGEN
37 template <int dim, int spacedim>
38 class Mapping;
39 
40 template <int dim>
41 class Quadrature;
42 
43 class ReferenceCell;
44 #endif
45 
46 
47 namespace internal
48 {
61  constexpr ReferenceCell
62  make_reference_cell_from_int(const std::uint8_t kind);
63 } // namespace internal
64 
65 
66 
98 {
99 public:
107  static ReferenceCell
108  n_vertices_to_type(const int dim, const unsigned int n_vertices);
109 
118  constexpr ReferenceCell();
119 
130  bool
131  is_hyper_cube() const;
132 
136  bool
137  is_simplex() const;
138 
143  unsigned int
144  get_dimension() const;
145 
159  template <int dim>
160  double
161  d_linear_shape_function(const Point<dim> &xi, const unsigned int i) const;
162 
167  template <int dim>
170  const unsigned int i) const;
171 
181  template <int dim, int spacedim = dim>
182  std::unique_ptr<Mapping<dim, spacedim>>
183  get_default_mapping(const unsigned int degree) const;
184 
196  template <int dim, int spacedim = dim>
197  const Mapping<dim, spacedim> &
199 
208  template <int dim>
210  get_gauss_type_quadrature(const unsigned n_points_1d) const;
211 
221  template <int dim>
223  get_midpoint_quadrature() const;
224 
238  template <int dim>
239  const Quadrature<dim> &
241 
256  unsigned int
257  n_vertices() const;
258 
264  vertex_indices() const;
265 
273  template <int dim>
274  Point<dim>
275  vertex(const unsigned int v) const;
276 
282  unsigned int
283  n_lines() const;
284 
290  line_indices() const;
291 
297  unsigned int
298  n_faces() const;
299 
305  face_indices() const;
306 
318  unsigned int
319  n_isotropic_children() const;
320 
326  isotropic_child_indices() const;
327 
341  face_reference_cell(const unsigned int face_no) const;
342 
356  static constexpr unsigned char
358 
369  static constexpr unsigned char
371 
413  unsigned int
414  child_cell_on_face(const unsigned int face,
415  const unsigned int subface,
416  const unsigned char face_orientation =
418 
428  std::array<unsigned int, 2>
429  standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const;
430 
440  std::array<unsigned int, 2>
441  standard_line_to_face_and_line_index(const unsigned int line) const;
442 
451  unsigned int
452  line_to_cell_vertices(const unsigned int line,
453  const unsigned int vertex) const;
454 
458  unsigned int
459  face_to_cell_lines(const unsigned int face,
460  const unsigned int line,
461  const unsigned char face_orientation) const;
462 
466  unsigned int
467  face_to_cell_vertices(const unsigned int face,
468  const unsigned int vertex,
469  const unsigned char face_orientation) const;
470 
490  template <int dim>
491  Point<dim>
492  face_vertex_location(const unsigned int face,
493  const unsigned int vertex) const;
494 
498  unsigned int
499  standard_to_real_face_vertex(const unsigned int vertex,
500  const unsigned int face,
501  const unsigned char face_orientation) const;
502 
506  unsigned int
507  standard_to_real_face_line(const unsigned int line,
508  const unsigned int face,
509  const unsigned char face_orientation) const;
510 
521  bool
522  standard_vs_true_line_orientation(const unsigned int line,
523  const unsigned char face_orientation,
524  const bool line_orientation) const;
525 
549  double
550  volume() const;
551 
564  template <int dim>
565  Point<dim>
566  barycenter() const;
567 
587  template <int dim>
588  bool
589  contains_point(const Point<dim> &p, const double tolerance = 0) const;
590 
595  template <int dim>
596  Point<dim>
597  closest_point(const Point<dim> &p) const;
598 
606  template <int dim>
608  unit_tangential_vectors(const unsigned int face_no,
609  const unsigned int i) const;
610 
614  template <int dim>
616  unit_normal_vectors(const unsigned int face_no) const;
617 
623  unsigned int
624  n_face_orientations(const unsigned int face_no) const;
625 
642  template <typename T, std::size_t N>
643  DEAL_II_DEPRECATED_EARLY unsigned char
644  compute_orientation(const std::array<T, N> &vertices_0,
645  const std::array<T, N> &vertices_1) const;
646 
672  template <typename T>
673  unsigned char
675  const ArrayView<const T> &vertices_1) const;
676 
677 
691  template <typename T, std::size_t N>
692  DEAL_II_DEPRECATED_EARLY std::array<T, N>
693  permute_according_orientation(const std::array<T, N> &vertices,
694  const unsigned int orientation) const;
695 
707  template <typename T>
708  boost::container::small_vector<T, 8>
710  const unsigned char orientation) const;
711 
716  faces_for_given_vertex(const unsigned int vertex_index) const;
717 
730  unsigned int
731  exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const;
732 
736  unsigned int
737  exodusii_face_to_deal_face(const unsigned int face_n) const;
738 
742  unsigned int
743  unv_vertex_to_deal_vertex(const unsigned int vertex_n) const;
744 
748  unsigned int
749  vtk_linear_type() const;
750 
755  unsigned int
756  vtk_quadratic_type() const;
757 
762  unsigned int
763  vtk_lagrange_type() const;
764 
777  template <int dim>
778  unsigned int
780  const std::array<unsigned, dim> &node_indices,
781  const std::array<unsigned, dim> &nodes_per_direction,
782  const bool legacy_format) const;
783 
787  unsigned int
788  vtk_vertex_to_deal_vertex(const unsigned int vertex_index) const;
789 
793  unsigned int
794  gmsh_element_type() const;
795 
809  std::string
810  to_string() const;
811 
815  constexpr operator std::uint8_t() const;
816 
820  constexpr bool
821  operator==(const ReferenceCell &type) const;
822 
826  constexpr bool
827  operator!=(const ReferenceCell &type) const;
828 
834  template <class Archive>
835  void
836  serialize(Archive &archive, const unsigned int /*version*/);
837 
841  static constexpr std::size_t
843 
848 private:
852  std::uint8_t kind;
853 
860  constexpr ReferenceCell(const std::uint8_t kind);
861 
866  {{{1, 0}}, {{0, 1}}}};
867 
868 
873  {{{0, 2, 1}},
874  {{0, 1, 2}},
875  {{2, 1, 0}},
876  {{1, 2, 0}},
877  {{1, 0, 2}},
878  {{2, 0, 1}}}};
879 
883  static constexpr ndarray<unsigned int, 8, 4>
885  {{0, 1, 2, 3}},
886  {{2, 3, 0, 1}},
887  {{2, 0, 3, 1}},
888  {{3, 1, 2, 0}},
889  {{3, 2, 1, 0}},
890  {{1, 0, 3, 2}},
891  {{1, 3, 0, 2}}}};
892 
897  friend constexpr ReferenceCell
899 
900  friend std::ostream &
901  operator<<(std::ostream &out, const ReferenceCell &reference_cell);
902 
903  friend std::istream &
904  operator>>(std::istream &in, ReferenceCell &reference_cell);
905 };
906 
907 
915 std::ostream &
916 operator<<(std::ostream &out, const ReferenceCell &reference_cell);
917 
925 std::istream &
926 operator>>(std::istream &in, ReferenceCell &reference_cell);
927 
928 
929 inline constexpr ReferenceCell::ReferenceCell(const std::uint8_t kind)
930  : kind(kind)
931 {}
932 
933 
934 
935 inline constexpr ReferenceCell::operator std::uint8_t() const
936 {
937  return kind;
938 }
939 
940 
941 
942 inline constexpr bool
944 {
945  return kind == type.kind;
946 }
947 
948 
949 
950 inline constexpr bool
952 {
953  return kind != type.kind;
954 }
955 
956 
957 
958 namespace internal
959 {
960  inline constexpr ReferenceCell
961  make_reference_cell_from_int(const std::uint8_t kind)
962  {
963 #ifndef DEAL_II_CXX14_CONSTEXPR_BUG
964  // Make sure these are the only indices from which objects can be
965  // created.
966  Assert((kind == std::numeric_limits<std::uint8_t>::max()) || (kind < 8),
967  ExcInternalError());
968 #endif
969 
970  // Call the private constructor, which we can from here because this
971  // function is a 'friend'.
972  return {kind};
973  }
974 } // namespace internal
975 
976 
977 
985 namespace ReferenceCells
986 {
987  constexpr const ReferenceCell Vertex =
989  constexpr const ReferenceCell Line =
991  constexpr const ReferenceCell Triangle =
993  constexpr const ReferenceCell Quadrilateral =
995  constexpr const ReferenceCell Tetrahedron =
997  constexpr const ReferenceCell Pyramid =
999  constexpr const ReferenceCell Wedge =
1001  constexpr const ReferenceCell Hexahedron =
1003  constexpr const ReferenceCell Invalid =
1006 
1012  template <int dim>
1013  constexpr const ReferenceCell &
1014  get_simplex();
1015 
1021  template <int dim>
1022  constexpr const ReferenceCell &
1023  get_hypercube();
1024 } // namespace ReferenceCells
1025 
1026 
1027 
1030 {}
1031 
1032 
1033 
1034 template <class Archive>
1035 inline void
1036 ReferenceCell::serialize(Archive &archive, const unsigned int /*version*/)
1037 {
1038  archive &kind;
1039 }
1040 
1041 
1042 
1043 inline constexpr std::size_t
1045 {
1046  return sizeof(ReferenceCells::Invalid);
1047 }
1048 
1049 
1050 
1052 ReferenceCell::faces_for_given_vertex(const unsigned int vertex) const
1053 {
1055  switch (this->kind)
1056  {
1057  case ReferenceCells::Line:
1058  return {&GeometryInfo<2>::vertex_to_face[vertex][0], 1};
1060  return {&GeometryInfo<2>::vertex_to_face[vertex][0], 2};
1062  {
1063  static constexpr ndarray<unsigned int, 3, 2> table = {
1064  {{{0, 2}}, {{0, 1}}, {{1, 2}}}};
1065  return table[vertex];
1066  }
1068  {
1069  static constexpr ndarray<unsigned int, 4, 3> table = {
1070  {{{0, 1, 2}}, {{0, 1, 3}}, {{0, 2, 3}}, {{1, 2, 3}}}};
1071 
1072  return table[vertex];
1073  }
1075  {
1076  static constexpr unsigned int X = numbers::invalid_unsigned_int;
1077  static constexpr ndarray<unsigned int, 5, 4> table = {
1078  {{{0, 1, 3, X}},
1079  {{0, 2, 3, X}},
1080  {{0, 1, 4, X}},
1081  {{0, 2, 4, X}},
1082  {{1, 2, 3, 4}}}};
1083 
1084  return {&table[vertex][0], vertex == 4 ? 4u : 3u};
1085  }
1086  case ReferenceCells::Wedge:
1087  {
1089  static constexpr ndarray<unsigned int, 6, 3> table = {{{{0, 2, 4}},
1090  {{0, 2, 3}},
1091  {{0, 3, 4}},
1092  {{1, 2, 4}},
1093  {{1, 2, 3}},
1094  {{1, 3, 4}}}};
1095 
1096  return table[vertex];
1097  }
1099  return {&GeometryInfo<3>::vertex_to_face[vertex][0], 3};
1100  default:
1101  Assert(false, ExcNotImplemented());
1102  }
1103 
1104  return {};
1105 }
1106 
1107 
1108 
1109 inline bool
1111 {
1112  return (*this == ReferenceCells::Vertex || *this == ReferenceCells::Line ||
1113  *this == ReferenceCells::Quadrilateral ||
1114  *this == ReferenceCells::Hexahedron);
1115 }
1116 
1117 
1118 
1119 inline bool
1121 {
1122  return (*this == ReferenceCells::Vertex || *this == ReferenceCells::Line ||
1123  *this == ReferenceCells::Triangle ||
1124  *this == ReferenceCells::Tetrahedron);
1125 }
1126 
1127 
1128 
1129 inline unsigned int
1131 {
1132  switch (this->kind)
1133  {
1135  return 0;
1136  case ReferenceCells::Line:
1137  return 1;
1140  return 2;
1143  case ReferenceCells::Wedge:
1145  return 3;
1146  default:
1147  Assert(false, ExcNotImplemented());
1148  }
1149 
1151 }
1152 
1153 
1154 
1155 template <int dim>
1158 {
1159  return Quadrature<dim>(std::vector<Point<dim>>({barycenter<dim>()}),
1160  std::vector<double>({volume()}));
1161 }
1162 
1163 
1164 
1165 inline unsigned int
1167 {
1168  switch (this->kind)
1169  {
1171  return 1;
1172  case ReferenceCells::Line:
1173  return 2;
1175  return 3;
1177  return 4;
1179  return 4;
1181  return 5;
1182  case ReferenceCells::Wedge:
1183  return 6;
1185  return 8;
1186  default:
1187  Assert(false, ExcNotImplemented());
1188  }
1189 
1191 }
1192 
1193 
1194 
1195 inline unsigned int
1197 {
1198  switch (this->kind)
1199  {
1201  return 0;
1202  case ReferenceCells::Line:
1203  return 1;
1205  return 3;
1207  return 4;
1209  return 6;
1211  return 8;
1212  case ReferenceCells::Wedge:
1213  return 9;
1215  return 12;
1216  default:
1217  Assert(false, ExcNotImplemented());
1218  }
1219 
1221 }
1222 
1223 
1224 
1225 template <int dim>
1226 Point<dim>
1227 ReferenceCell::vertex(const unsigned int v) const
1228 {
1231 
1232  switch (dim)
1233  {
1234  case 0:
1235  {
1236  if (*this == ReferenceCells::Vertex)
1237  return Point<dim>();
1238  break;
1239  }
1240  case 1:
1241  {
1242  static const Point<dim> vertices[2] = {
1243  Point<dim>(), // the origin
1244  Point<dim>::unit_vector(0) // unit point along x-axis
1245  };
1246  if (*this == ReferenceCells::Line)
1247  return vertices[v];
1248  break;
1249  }
1250  case 2:
1251  {
1252  switch (this->kind)
1253  {
1255  {
1256  static const Point<dim> vertices[3] = {
1257  Point<dim>(), // the origin
1258  Point<dim>::unit_vector(0), // unit point along x-axis
1259  Point<dim>::unit_vector(1) // unit point along y-axis
1260  };
1261  return vertices[v];
1262  }
1264  {
1265  static const Point<dim> vertices[4] = {
1266  // First the two points on the x-axis
1267  Point<dim>(),
1269  // Then these two points shifted in the y-direction
1272  return vertices[v];
1273  }
1274  }
1275  break;
1276  }
1277  case 3:
1278  {
1279  switch (this->kind)
1280  {
1282  {
1283  static const Point<dim> vertices[4] = {
1284  Point<dim>(), // the origin
1285  Point<dim>::unit_vector(0), // unit point along x-axis
1286  Point<dim>::unit_vector(1), // unit point along y-axis
1287  Point<dim>::unit_vector(2) // unit point along z-axis
1288  };
1289  return vertices[v];
1290  }
1292  {
1293  static const Point<dim> vertices[5] = {
1294  Point<dim>{-1.0, -1.0, 0.0},
1295  Point<dim>{+1.0, -1.0, 0.0},
1296  Point<dim>{-1.0, +1.0, 0.0},
1297  Point<dim>{+1.0, +1.0, 0.0},
1298  Point<dim>{+0.0, +0.0, 1.0}};
1299  return vertices[v];
1300  }
1301  case ReferenceCells::Wedge:
1302  {
1303  static const Point<dim> vertices[6] = {
1304  // First the three points on the triangular base of the
1305  // wedge:
1306  Point<dim>(),
1309  // And now everything shifted in the z-direction again
1313  return vertices[v];
1314  }
1316  {
1317  static const Point<dim> vertices[8] = {
1318  // First the two points on the x-axis
1319  Point<dim>(),
1321  // Then these two points shifted in the y-direction
1324  // And now all four points shifted in the z-direction
1331  return vertices[v];
1332  }
1333  }
1334  break;
1335  }
1336  default:
1337  Assert(false, ExcNotImplemented());
1338  }
1339 
1340  Assert(false, ExcNotImplemented());
1341  return Point<dim>();
1342 }
1343 
1344 
1345 inline unsigned int
1347 {
1348  switch (this->kind)
1349  {
1351  return 0;
1352  case ReferenceCells::Line:
1353  return 2;
1355  return 3;
1357  return 4;
1359  return 4;
1361  return 5;
1362  case ReferenceCells::Wedge:
1363  return 5;
1365  return 6;
1366  default:
1367  Assert(false, ExcNotImplemented());
1368  }
1369 
1371 }
1372 
1373 
1374 
1377 {
1378  return {0U, n_faces()};
1379 }
1380 
1381 
1382 
1383 inline unsigned int
1385 {
1386  switch (this->kind)
1387  {
1389  return 0;
1390  case ReferenceCells::Line:
1391  return 2;
1393  return 4;
1395  return 4;
1397  return 8;
1399  // We haven't yet decided how to refine pyramids. Update this when we
1400  // have
1401  Assert(false, ExcNotImplemented());
1403  case ReferenceCells::Wedge:
1404  return 8;
1406  return 8;
1407  default:
1408  Assert(false, ExcNotImplemented());
1409  }
1410 
1412 }
1413 
1414 
1415 
1418 {
1419  return {0U, n_isotropic_children()};
1420 }
1421 
1422 
1423 
1426 {
1427  return {0U, n_vertices()};
1428 }
1429 
1430 
1431 
1434 {
1435  return {0U, n_lines()};
1436 }
1437 
1438 
1439 
1440 inline ReferenceCell
1441 ReferenceCell::face_reference_cell(const unsigned int face_no) const
1442 {
1443  AssertIndexRange(face_no, n_faces());
1444 
1445  switch (this->kind)
1446  {
1448  return ReferenceCells::Invalid;
1449  case ReferenceCells::Line:
1450  return ReferenceCells::Vertex;
1453  return ReferenceCells::Line;
1455  return ReferenceCells::Triangle;
1457  if (face_no == 0)
1459  else
1460  return ReferenceCells::Triangle;
1461  case ReferenceCells::Wedge:
1462  if (face_no > 1)
1464  else
1465  return ReferenceCells::Triangle;
1468  default:
1469  Assert(false, ExcNotImplemented());
1470  }
1471 
1472  return ReferenceCells::Invalid;
1473 }
1474 
1475 
1476 
1477 inline constexpr unsigned char
1479 {
1480  // Our convention is that 'orientation' has a default value of true and
1481  // occupies the least-significant bit while rotate and flip have default
1482  // values of 'false' and occupy the second and third bits.
1483  return 0b001;
1484 }
1485 
1486 
1487 
1488 inline constexpr unsigned char
1490 {
1491  // For a reversed line 'orientation' is false and neither flip nor rotate are
1492  // defined.
1493  return 0b000;
1494 }
1495 
1496 
1497 
1498 inline unsigned int
1500  const unsigned int face,
1501  const unsigned int subface,
1502  const unsigned char combined_face_orientation) const
1503 {
1504  AssertIndexRange(face, n_faces());
1506 
1507  switch (this->kind)
1508  {
1510  case ReferenceCells::Line:
1511  {
1512  Assert(false, ExcNotImplemented());
1513  break;
1514  }
1516  {
1517  static constexpr ndarray<unsigned int, 3, 2> subcells = {
1518  {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
1519 
1520  return subcells[face][subface];
1521  }
1523  {
1524  const bool face_orientation =
1525  Utilities::get_bit(combined_face_orientation, 0);
1526  const bool face_flip =
1527  Utilities::get_bit(combined_face_orientation, 2);
1528  const bool face_rotation =
1529  Utilities::get_bit(combined_face_orientation, 1);
1530 
1533  face,
1534  subface,
1535  face_orientation,
1536  face_flip,
1537  face_rotation);
1538  }
1541  case ReferenceCells::Wedge:
1542  {
1543  Assert(false, ExcNotImplemented());
1544  break;
1545  }
1547  {
1548  const bool face_orientation =
1549  Utilities::get_bit(combined_face_orientation, 0);
1550  const bool face_flip =
1551  Utilities::get_bit(combined_face_orientation, 2);
1552  const bool face_rotation =
1553  Utilities::get_bit(combined_face_orientation, 1);
1554 
1557  face,
1558  subface,
1559  face_orientation,
1560  face_flip,
1561  face_rotation);
1562  }
1563  default:
1564  Assert(false, ExcNotImplemented());
1565  }
1566 
1568 }
1569 
1570 
1571 
1572 inline std::array<unsigned int, 2>
1574  const unsigned int vertex) const
1575 {
1577  // Work around a GCC warning at higher optimization levels by making all of
1578  // these tables the same size
1579  constexpr unsigned int X = numbers::invalid_unsigned_int;
1580 
1581  switch (this->kind)
1582  {
1584  case ReferenceCells::Line:
1585  Assert(false, ExcNotImplemented());
1586  break;
1588  {
1589  static constexpr ndarray<unsigned int, 6, 2> table = {
1590  {{{0, 0}}, {{0, 1}}, {{1, 1}}, {{X, X}}, {{X, X}}, {{X, X}}}};
1591 
1592  return table[vertex];
1593  }
1595  {
1597  vertex);
1598  }
1600  {
1601  static constexpr ndarray<unsigned int, 6, 2> table = {
1602  {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 2}}, {{X, X}}, {{X, X}}}};
1603 
1604  return table[vertex];
1605  }
1607  {
1608  static constexpr ndarray<unsigned int, 6, 2> table = {
1609  {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{0, 3}}, {{1, 2}}, {{X, X}}}};
1610 
1611  return table[vertex];
1612  }
1613  case ReferenceCells::Wedge:
1614  {
1615  static constexpr ndarray<unsigned int, 6, 2> table = {
1616  {{{0, 1}}, {{0, 0}}, {{0, 2}}, {{1, 0}}, {{1, 1}}, {{1, 2}}}};
1617 
1618  return table[vertex];
1619  }
1621  {
1623  vertex);
1624  }
1625  default:
1626  Assert(false, ExcNotImplemented());
1627  }
1628 
1629  return {};
1630 }
1631 
1632 
1633 
1634 inline std::array<unsigned int, 2>
1636  const unsigned int line) const
1637 {
1638  AssertIndexRange(line, n_lines());
1639 
1640  switch (this->kind)
1641  {
1643  case ReferenceCells::Line:
1646  {
1647  Assert(false, ExcNotImplemented());
1648  break;
1649  }
1651  {
1652  static const std::array<unsigned int, 2> table[6] = {
1653  {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}};
1654 
1655  return table[line];
1656  }
1658  {
1659  static const std::array<unsigned int, 2> table[8] = {{{0, 0}},
1660  {{0, 1}},
1661  {{0, 2}},
1662  {{0, 3}},
1663  {{1, 2}},
1664  {{2, 1}},
1665  {{1, 1}},
1666  {{2, 2}}};
1667 
1668  return table[line];
1669  }
1670  case ReferenceCells::Wedge:
1671  {
1672  static const std::array<unsigned int, 2> table[9] = {{{0, 0}},
1673  {{0, 2}},
1674  {{0, 1}},
1675  {{1, 0}},
1676  {{1, 1}},
1677  {{1, 2}},
1678  {{2, 0}},
1679  {{2, 1}},
1680  {{3, 1}}};
1681 
1682  return table[line];
1683  }
1685  {
1687  }
1688  default:
1689  Assert(false, ExcNotImplemented());
1690  }
1691 
1692  return {};
1693 }
1694 
1695 
1696 inline unsigned int
1698  const unsigned int vertex) const
1699 {
1701  AssertIndexRange(line, n_lines());
1702 
1703  switch (this->kind)
1704  {
1706  case ReferenceCells::Line:
1707  return vertex;
1709  {
1710  static constexpr ndarray<unsigned int, 3, 2> table = {
1711  {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
1712  return table[line][vertex];
1713  }
1715  {
1716  static constexpr ndarray<unsigned int, 4, 2> table = {
1717  {{{0, 2}}, {{1, 3}}, {{0, 1}}, {{2, 3}}}};
1718  return table[line][vertex];
1719  }
1721  {
1722  static constexpr ndarray<unsigned int, 6, 2> table = {
1723  {{{0, 1}}, {{1, 2}}, {{2, 0}}, {{0, 3}}, {{1, 3}}, {{2, 3}}}};
1724  return table[line][vertex];
1725  }
1727  {
1728  static constexpr ndarray<unsigned int, 8, 2> table = {{{{0, 2}},
1729  {{1, 3}},
1730  {{0, 1}},
1731  {{2, 3}},
1732  {{4, 0}},
1733  {{1, 4}},
1734  {{2, 4}},
1735  {{4, 3}}}};
1736  return table[line][vertex];
1737  }
1738  case ReferenceCells::Wedge:
1739  {
1740  static constexpr ndarray<unsigned int, 9, 2> table = {{{{1, 0}},
1741  {{2, 1}},
1742  {{0, 2}},
1743  {{3, 4}},
1744  {{4, 5}},
1745  {{5, 3}},
1746  {{0, 3}},
1747  {{1, 4}},
1748  {{2, 5}}}};
1749  return table[line][vertex];
1750  }
1752  {
1753  // first four lines comprise the bottom face, next four are the top,
1754  // and the last four are 'bottom to top'
1755  static constexpr ndarray<unsigned int, 12, 2> table = {{{{0, 2}},
1756  {{1, 3}},
1757  {{0, 1}},
1758  {{2, 3}},
1759  {{4, 6}},
1760  {{5, 7}},
1761  {{4, 5}},
1762  {{6, 7}},
1763  {{0, 4}},
1764  {{1, 5}},
1765  {{2, 6}},
1766  {{3, 7}}}};
1767  return table[line][vertex];
1768  }
1769 
1770  default:
1771  Assert(false, ExcNotImplemented());
1772  }
1773 
1775 }
1776 
1777 
1778 inline unsigned int
1779 ReferenceCell::face_to_cell_lines(const unsigned int face,
1780  const unsigned int line,
1781  const unsigned char face_orientation) const
1782 {
1783  AssertIndexRange(face, n_faces());
1785 
1786  static constexpr unsigned int X = numbers::invalid_unsigned_int;
1787 
1788  switch (this->kind)
1789  {
1791  {
1792  Assert(false, ExcNotImplemented());
1793  break;
1794  }
1795  case ReferenceCells::Line:
1796  {
1798  face,
1799  line,
1800  Utilities::get_bit(face_orientation, 0),
1801  Utilities::get_bit(face_orientation, 2),
1802  Utilities::get_bit(face_orientation, 1));
1803  }
1805  {
1806  return face;
1807  }
1809  {
1811  face,
1812  line,
1813  Utilities::get_bit(face_orientation, 0),
1814  Utilities::get_bit(face_orientation, 2),
1815  Utilities::get_bit(face_orientation, 1));
1816  }
1818  {
1819  static constexpr ndarray<unsigned int, 4, 3> table = {
1820  {{{0, 1, 2}}, {{0, 3, 4}}, {{2, 5, 3}}, {{1, 4, 5}}}};
1821 
1822  return table[face][standard_to_real_face_line(
1823  line, face, face_orientation)];
1824  }
1826  {
1827  static constexpr ndarray<unsigned int, 5, 4> table = {
1828  {{{0, 1, 2, 3}},
1829  {{0, 6, 4, X}},
1830  {{1, 5, 7, X}},
1831  {{2, 4, 5, X}},
1832  {{3, 7, 6, X}}}};
1833 
1834  return table[face][standard_to_real_face_line(
1835  line, face, face_orientation)];
1836  }
1837  case ReferenceCells::Wedge:
1838  {
1839  static constexpr ndarray<unsigned int, 5, 4> table = {
1840  {{{0, 2, 1, X}},
1841  {{3, 4, 5, X}},
1842  {{6, 7, 0, 3}},
1843  {{7, 8, 1, 4}},
1844  {{8, 6, 5, 2}}}};
1845 
1846  return table[face][standard_to_real_face_line(
1847  line, face, face_orientation)];
1848  }
1850  {
1852  face,
1853  line,
1854  Utilities::get_bit(face_orientation, 0),
1855  Utilities::get_bit(face_orientation, 2),
1856  Utilities::get_bit(face_orientation, 1));
1857  }
1858  default:
1859  Assert(false, ExcNotImplemented());
1860  }
1861 
1863 }
1864 
1865 
1866 
1867 inline unsigned int
1869  const unsigned int vertex,
1870  const unsigned char face_orientation) const
1871 {
1872  AssertIndexRange(face, n_faces());
1874 
1875  switch (this->kind)
1876  {
1878  {
1879  Assert(false, ExcNotImplemented());
1880  break;
1881  }
1882  case ReferenceCells::Line:
1883  {
1885  face,
1886  vertex,
1887  Utilities::get_bit(face_orientation, 0),
1888  Utilities::get_bit(face_orientation, 2),
1889  Utilities::get_bit(face_orientation, 1));
1890  }
1892  {
1893  static constexpr ndarray<unsigned int, 3, 2> table = {
1894  {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
1895 
1896  return table[face][face_orientation != 0u ? vertex : (1 - vertex)];
1897  }
1899  {
1901  face,
1902  vertex,
1903  Utilities::get_bit(face_orientation, 0),
1904  Utilities::get_bit(face_orientation, 2),
1905  Utilities::get_bit(face_orientation, 1));
1906  }
1908  {
1909  static constexpr ndarray<unsigned int, 4, 3> table = {
1910  {{{0, 1, 2}}, {{1, 0, 3}}, {{0, 2, 3}}, {{2, 1, 3}}}};
1911 
1912  return table[face][standard_to_real_face_vertex(
1913  vertex, face, face_orientation)];
1914  }
1916  {
1917  constexpr auto X = numbers::invalid_unsigned_int;
1918  static constexpr ndarray<unsigned int, 5, 4> table = {
1919  {{{0, 1, 2, 3}},
1920  {{0, 2, 4, X}},
1921  {{3, 1, 4, X}},
1922  {{1, 0, 4, X}},
1923  {{2, 3, 4, X}}}};
1924 
1925  return table[face][standard_to_real_face_vertex(
1926  vertex, face, face_orientation)];
1927  }
1928  case ReferenceCells::Wedge:
1929  {
1930  constexpr auto X = numbers::invalid_unsigned_int;
1931  static constexpr ndarray<unsigned int, 6, 4> table = {
1932  {{{1, 0, 2, X}},
1933  {{3, 4, 5, X}},
1934  {{0, 1, 3, 4}},
1935  {{1, 2, 4, 5}},
1936  {{2, 0, 5, 3}}}};
1937 
1938  return table[face][standard_to_real_face_vertex(
1939  vertex, face, face_orientation)];
1940  }
1942  {
1944  face,
1945  vertex,
1946  Utilities::get_bit(face_orientation, 0),
1947  Utilities::get_bit(face_orientation, 2),
1948  Utilities::get_bit(face_orientation, 1));
1949  }
1950  default:
1951  Assert(false, ExcNotImplemented());
1952  }
1953 
1955 }
1956 
1957 
1958 
1959 template <int dim>
1960 Point<dim>
1961 ReferenceCell::face_vertex_location(const unsigned int face,
1962  const unsigned int vertex) const
1963 {
1964  AssertIndexRange(face, n_faces());
1966  AssertDimension(dim, this->get_dimension());
1967 
1968  switch (this->kind)
1969  {
1971  {
1972  Assert(false, ExcNotImplemented());
1973  break;
1974  }
1975  case ReferenceCells::Line:
1976  {
1977  if (face == 0)
1978  return Point<dim>(0);
1979  else
1980  return Point<dim>(1);
1981  }
1983  {
1984  static const ndarray<Point<dim>, 3, 2> table = {
1985  {{{Point<dim>(0, 0), Point<dim>(1, 0)}},
1986  {{Point<dim>(1, 0), Point<dim>(0, 1)}},
1987  {{Point<dim>(0, 1), Point<dim>(0, 0)}}}};
1988 
1989  return table[face][vertex];
1990  }
1992  {
1993  static const ndarray<Point<dim>, 4, 2> table = {
1994  {{{Point<dim>(0, 0), Point<dim>(0, 1)}},
1995  {{Point<dim>(1, 0), Point<dim>(1, 1)}},
1996  {{Point<dim>(0, 0), Point<dim>(1, 0)}},
1997  {{Point<dim>(0, 1), Point<dim>(1, 1)}}}};
1998 
1999  return table[face][vertex];
2000  }
2002  {
2003  static const ndarray<Point<dim>, 4, 3> table = {
2004  {{{Point<dim>(0.0, 0.0, 0.0),
2005  Point<dim>(1.0, 0.0, 0.0),
2006  Point<dim>(0.0, 1.0, 0.0)}},
2007  {{Point<dim>(1.0, 0.0, 0.0),
2008  Point<dim>(0.0, 0.0, 0.0),
2009  Point<dim>(0.0, 0.0, 1.0)}},
2010  {{Point<dim>(0.0, 0.0, 0.0),
2011  Point<dim>(0.0, 1.0, 0.0),
2012  Point<dim>(0.0, 0.0, 1.0)}},
2013  {{Point<dim>(0.0, 1.0, 0.0),
2014  Point<dim>(1.0, 0.0, 0.0),
2015  Point<dim>(0.0, 0.0, 1.0)}}}};
2016 
2017  return table[face][vertex];
2018  }
2020  {
2021  // We would like to use ndarray here as we do above, but
2022  // not every face has the same number of vertices and so
2023  // the inner array needs to be a dynamic vector:
2024  static const std::array<std::vector<Point<dim>>, 5> table = {
2025  {{{Point<dim>(-1.0, -1.0, 0.0),
2026  Point<dim>(+1.0, -1.0, 0.0),
2027  Point<dim>(-1.0, +1.0, 0.0),
2028  Point<dim>(+1.0, +1.0, 0.0)}},
2029  {{Point<dim>(-1.0, -1.0, 0.0),
2030  Point<dim>(-1.0, +1.0, 0.0),
2031  Point<dim>(+0.0, +0.0, 1.0)}},
2032  {{Point<dim>(+1.0, +1.0, 0.0),
2033  Point<dim>(+1.0, -1.0, 0.0),
2034  Point<dim>(+0.0, +0.0, 1.0)}},
2035  {{Point<dim>(+1.0, -1.0, 0.0),
2036  Point<dim>(-1.0, -1.0, 0.0),
2037  Point<dim>(+0.0, +0.0, 1.0)}},
2038  {{Point<dim>(-1.0, +1.0, 0.0),
2039  Point<dim>(+1.0, +1.0, 0.0),
2040  Point<dim>(+0.0, +0.0, 1.0)}}}};
2041 
2042  return table[face][vertex];
2043  }
2044  case ReferenceCells::Wedge:
2045  {
2046  // We would like to use ndarray here as we do above, but
2047  // not every face has the same number of vertices and so
2048  // the inner array needs to be a dynamic vector:
2049  static const std::array<std::vector<Point<dim>>, 5> table = {
2050  {{{Point<dim>(1.0, 0.0, 0.0),
2051  Point<dim>(0.0, 0.0, 0.0),
2052  Point<dim>(0.0, 1.0, 0.0)}},
2053  {{Point<dim>(0.0, 0.0, 1.0),
2054  Point<dim>(1.0, 0.0, 1.0),
2055  Point<dim>(0.0, 1.0, 1.0)}},
2056  {{Point<dim>(0.0, 0.0, 0.0),
2057  Point<dim>(1.0, 0.0, 0.0),
2058  Point<dim>(0.0, 0.0, 1.0),
2059  Point<dim>(1.0, 0.0, 1.0)}},
2060  {{Point<dim>(1.0, 0.0, 0.0),
2061  Point<dim>(0.0, 1.0, 0.0),
2062  Point<dim>(1.0, 0.0, 1.0),
2063  Point<dim>(0.0, 1.0, 1.0)}},
2064  {{Point<dim>(0.0, 1.0, 0.0),
2065  Point<dim>(0.0, 0.0, 0.0),
2066  Point<dim>(0.0, 1.0, 1.0),
2067  Point<dim>(0.0, 0.0, 1.0)}}}};
2068 
2069  return table[face][vertex];
2070  }
2072  {
2073  /*
2074  * Follow these figures:
2075  *
2076  * * <li> Faces 0 and 1:
2077  * @verbatim
2078  * Face 0 Face 1
2079  * *-------* *-------*
2080  * /| | / /|
2081  * 3 1 | / 3 1
2082  * y/ | | / y/ |
2083  * * |x | *-------* |x
2084  * | *-------* | | *
2085  * 0 / / | 0 /
2086  * | 2 / | | 2
2087  * |/ / | |/
2088  * *-------* *-------*
2089  * @endverbatim
2090  *
2091  * <li> Faces 2 and 3:
2092  * @verbatim
2093  * x Face 3 Face 2
2094  * *---1---* *-------*
2095  * /| | / /|
2096  * / | 3 / / |
2097  * / 2 | x/ / |
2098  * * | | *---1---* |
2099  * | *---0---*y | | *
2100  * | / / | 3 /
2101  * | / / 2 | /
2102  * |/ / | |/
2103  * *-------* *---0---*y
2104  * @endverbatim
2105  *
2106  * <li> Faces 4 and 5:
2107  * @verbatim
2108  * Face 4 y Face 5
2109  * *-------* *---3---*
2110  * /| | / /|
2111  * / | | 0 1 |
2112  * / | | / / |
2113  * * |y | *---2---* x |
2114  * | *---3---* | | *
2115  * | / / | | /
2116  * | 0 1 | | /
2117  * |/ / | |/
2118  * *---2---* x *-------*
2119  * @endverbatim
2120  * </ul>
2121  */
2122  static const ndarray<Point<dim>, 6, 4> table = {
2123  {// Face 0:
2124  {{Point<dim>(0, 0, 0),
2125  Point<dim>(0, 1, 0),
2126  Point<dim>(0, 0, 1),
2127  Point<dim>(0, 1, 1)}},
2128  // Face 1:
2129  {{Point<dim>(1, 0, 0),
2130  Point<dim>(1, 1, 0),
2131  Point<dim>(1, 0, 1),
2132  Point<dim>(1, 1, 1)}},
2133  // Face 2:
2134  {{Point<dim>(0, 0, 0),
2135  Point<dim>(0, 0, 1),
2136  Point<dim>(1, 0, 0),
2137  Point<dim>(1, 0, 1)}},
2138  // Face 3:
2139  {{Point<dim>(0, 1, 0),
2140  Point<dim>(0, 1, 1),
2141  Point<dim>(1, 1, 0),
2142  Point<dim>(1, 1, 1)}},
2143  // Face 4:
2144  {{Point<dim>(0, 0, 0),
2145  Point<dim>(1, 0, 0),
2146  Point<dim>(0, 1, 0),
2147  Point<dim>(1, 1, 0)}},
2148  // Face 5:
2149  {{Point<dim>(0, 0, 1),
2150  Point<dim>(1, 0, 1),
2151  Point<dim>(0, 1, 1),
2152  Point<dim>(1, 1, 1)}}}};
2153 
2154  return table[face][vertex];
2155  }
2156  default:
2157  Assert(false, ExcNotImplemented());
2158  }
2159 
2160  return {};
2161 }
2162 
2163 
2164 
2165 inline unsigned int
2167  const unsigned int vertex,
2168  const unsigned int face,
2169  const unsigned char face_orientation) const
2170 {
2171  AssertIndexRange(face, n_faces());
2173 
2174  switch (this->kind)
2175  {
2177  case ReferenceCells::Line:
2178  Assert(false, ExcNotImplemented());
2179  break;
2182  return line_vertex_permutations[face_orientation][vertex];
2184  return triangle_vertex_permutations[face_orientation][vertex];
2186  // face 0 is a quadrilateral
2187  if (face == 0)
2188  return quadrilateral_vertex_permutations[face_orientation][vertex];
2189  else
2190  return triangle_vertex_permutations[face_orientation][vertex];
2191  case ReferenceCells::Wedge:
2192  // faces 0 and 1 are triangles
2193  if (face > 1)
2194  return quadrilateral_vertex_permutations[face_orientation][vertex];
2195  else
2196  return triangle_vertex_permutations[face_orientation][vertex];
2198  return quadrilateral_vertex_permutations[face_orientation][vertex];
2199  default:
2200  Assert(false, ExcNotImplemented());
2201  }
2202 
2203  Assert(false, ExcNotImplemented());
2205 }
2206 
2207 
2208 
2209 inline unsigned int
2211  const unsigned int line,
2212  const unsigned int face,
2213  const unsigned char face_orientation) const
2214 {
2215  AssertIndexRange(face, n_faces());
2217 
2218  static constexpr ndarray<unsigned int, 6, 3> triangle_table = {{{{2, 1, 0}},
2219  {{0, 1, 2}},
2220  {{1, 0, 2}},
2221  {{1, 2, 0}},
2222  {{0, 2, 1}},
2223  {{2, 0, 1}}}};
2224 
2225 
2226  switch (this->kind)
2227  {
2229  case ReferenceCells::Line:
2232  Assert(false, ExcNotImplemented());
2233  break;
2235  return triangle_table[face_orientation][line];
2237  if (face == 0) // The quadrilateral face
2238  {
2240  line,
2241  Utilities::get_bit(face_orientation, 0),
2242  Utilities::get_bit(face_orientation, 2),
2243  Utilities::get_bit(face_orientation, 1));
2244  }
2245  else // One of the triangular faces
2246  {
2247  return triangle_table[face_orientation][line];
2248  }
2249  case ReferenceCells::Wedge:
2250  if (face > 1) // One of the quadrilateral faces
2251  {
2253  line,
2254  Utilities::get_bit(face_orientation, 0),
2255  Utilities::get_bit(face_orientation, 2),
2256  Utilities::get_bit(face_orientation, 1));
2257  }
2258  else // One of the triangular faces
2259  return triangle_table[face_orientation][line];
2261  {
2262  static constexpr ndarray<unsigned int, 8, 4> table = {
2263  {{{2, 3, 0, 1}},
2264  {{0, 1, 2, 3}},
2265  {{0, 1, 3, 2}},
2266  {{3, 2, 0, 1}},
2267  {{3, 2, 1, 0}},
2268  {{1, 0, 3, 2}},
2269  {{1, 0, 2, 3}},
2270  {{2, 3, 1, 0}}}};
2271  return table[face_orientation][line];
2272  }
2273  default:
2274  Assert(false, ExcNotImplemented());
2275  }
2276 
2278 }
2279 
2280 
2281 
2282 namespace ReferenceCells
2283 {
2284  template <int dim>
2285  inline constexpr const ReferenceCell &
2287  {
2288  switch (dim)
2289  {
2290  case 0:
2291  return ReferenceCells::Vertex;
2292  case 1:
2293  return ReferenceCells::Line;
2294  case 2:
2295  return ReferenceCells::Triangle;
2296  case 3:
2298  default:
2299  Assert(false, ExcNotImplemented());
2300  }
2301  return ReferenceCells::Invalid;
2302  }
2303 
2304 
2305 
2306  template <int dim>
2307  inline constexpr const ReferenceCell &
2309  {
2310  switch (dim)
2311  {
2312  case 0:
2313  return ReferenceCells::Vertex;
2314  case 1:
2315  return ReferenceCells::Line;
2316  case 2:
2318  case 3:
2320  default:
2321  Assert(false, ExcNotImplemented());
2322  }
2323  return ReferenceCells::Invalid;
2324  }
2325 } // namespace ReferenceCells
2326 
2327 
2328 inline ReferenceCell
2329 ReferenceCell::n_vertices_to_type(const int dim, const unsigned int n_vertices)
2330 {
2331  AssertIndexRange(dim, 4);
2333 
2334  const auto X = ReferenceCells::Invalid;
2335  static const ndarray<ReferenceCell, 4, 9> table = {
2336  {// dim 0
2337  {{X, ReferenceCells::Vertex, X, X, X, X, X, X, X}},
2338  // dim 1
2339  {{X, X, ReferenceCells::Line, X, X, X, X, X, X}},
2340  // dim 2
2341  {{X,
2342  X,
2343  X,
2346  X,
2347  X,
2348  X,
2349  X}},
2350  // dim 3
2351  {{X,
2352  X,
2353  X,
2354  X,
2358  X,
2360  Assert(table[dim][n_vertices] != ReferenceCells::Invalid,
2361  ExcMessage("The combination of dim = " + std::to_string(dim) +
2362  " and n_vertices = " + std::to_string(n_vertices) +
2363  " does not correspond to a known reference cell type."));
2364  return table[dim][n_vertices];
2365 }
2366 
2367 
2368 
2369 template <int dim>
2370 inline double
2372  const unsigned int i) const
2373 {
2376  switch (this->kind)
2377  {
2379  case ReferenceCells::Line:
2383  // see also BarycentricPolynomials<2>::compute_value
2385  {
2386  switch (i)
2387  {
2388  case 0:
2389  return 1.0 - xi[std::min(0, dim - 1)] -
2390  xi[std::min(1, dim - 1)];
2391  case 1:
2392  return xi[std::min(0, dim - 1)];
2393  case 2:
2394  return xi[std::min(1, dim - 1)];
2395  default:
2396  Assert(false, ExcInternalError());
2397  }
2398  }
2399  // see also BarycentricPolynomials<3>::compute_value
2401  {
2402  switch (i)
2403  {
2404  case 0:
2405  return 1.0 - xi[std::min(0, dim - 1)] -
2406  xi[std::min(1, dim - 1)] - xi[std::min(2, dim - 1)];
2407  case 1:
2408  return xi[std::min(0, dim - 1)];
2409  case 2:
2410  return xi[std::min(1, dim - 1)];
2411  case 3:
2412  return xi[std::min(2, dim - 1)];
2413  default:
2414  Assert(false, ExcInternalError());
2415  }
2416  }
2417  // see also ScalarLagrangePolynomialPyramid::compute_value()
2419  {
2420  const double Q14 = 0.25;
2421 
2422  const double r = xi[std::min(0, dim - 1)];
2423  const double s = xi[std::min(1, dim - 1)];
2424  const double t = xi[std::min(2, dim - 1)];
2425 
2426  const double ratio =
2427  (std::fabs(t - 1.0) > 1.0e-14 ? (r * s * t) / (1.0 - t) : 0.0);
2428 
2429  if (i == 0)
2430  return Q14 * ((1.0 - r) * (1.0 - s) - t + ratio);
2431  if (i == 1)
2432  return Q14 * ((1.0 + r) * (1.0 - s) - t - ratio);
2433  if (i == 2)
2434  return Q14 * ((1.0 - r) * (1.0 + s) - t - ratio);
2435  if (i == 3)
2436  return Q14 * ((1.0 + r) * (1.0 + s) - t + ratio);
2437  else
2438  return t;
2439  }
2440  // see also ScalarLagrangePolynomialWedge::compute_value()
2441  case ReferenceCells::Wedge:
2443  .d_linear_shape_function<2>(Point<2>(xi[std::min(0, dim - 1)],
2444  xi[std::min(1, dim - 1)]),
2445  i % 3) *
2447  .d_linear_shape_function<1>(Point<1>(xi[std::min(2, dim - 1)]),
2448  i / 3);
2449  default:
2450  Assert(false, ExcNotImplemented());
2451  }
2452 
2453  return 0.0;
2454 }
2455 
2456 
2457 
2458 template <int dim>
2459 inline Tensor<1, dim>
2461  const unsigned int i) const
2462 {
2464  switch (this->kind)
2465  {
2467  case ReferenceCells::Line:
2471  // see also BarycentricPolynomials<2>::compute_grad()
2473  switch (i)
2474  {
2475  case 0:
2476  return Point<dim>(-1.0, -1.0);
2477  case 1:
2478  return Point<dim>(+1.0, +0.0);
2479  case 2:
2480  return Point<dim>(+0.0, +1.0);
2481  default:
2482  Assert(false, ExcInternalError());
2483  }
2484  default:
2485  Assert(false, ExcNotImplemented());
2486  }
2487 
2488  return Point<dim>(+0.0, +0.0, +0.0);
2489 }
2490 
2491 
2492 
2493 inline double
2495 {
2496  switch (this->kind)
2497  {
2499  return 0;
2500  case ReferenceCells::Line:
2501  return 1;
2503  return 1. / 2.;
2505  return 1;
2507  return 1. / 6.;
2509  return 4. / 3.;
2510  case ReferenceCells::Wedge:
2511  return 1. / 2.;
2513  return 1;
2514  default:
2515  Assert(false, ExcNotImplemented());
2516  }
2517 
2518  return 0.0;
2519 }
2520 
2521 
2522 
2523 template <int dim>
2524 inline Point<dim>
2526 {
2528 
2529  switch (this->kind)
2530  {
2532  return Point<dim>();
2533  case ReferenceCells::Line:
2534  return Point<dim>(1. / 2.);
2536  return Point<dim>(1. / 3., 1. / 3.);
2538  return Point<dim>(1. / 2., 1. / 2.);
2540  return Point<dim>(1. / 4., 1. / 4., 1. / 4.);
2542  return Point<dim>(0, 0, 1. / 4.);
2543  case ReferenceCells::Wedge:
2544  return Point<dim>(1. / 3, 1. / 3, 1. / 2.);
2546  return Point<dim>(1. / 2., 1. / 2., 1. / 2.);
2547  default:
2548  Assert(false, ExcNotImplemented());
2549  }
2550 
2551  return Point<dim>();
2552 }
2553 
2554 
2555 
2556 template <int dim>
2557 inline bool
2558 ReferenceCell::contains_point(const Point<dim> &p, const double tolerance) const
2559 {
2561 
2562  // Introduce abbreviations to silence compiler warnings about invalid
2563  // array accesses (that can't happen, of course, but the compiler
2564  // doesn't know that).
2565  constexpr unsigned int x_coordinate = 0;
2566  constexpr unsigned int y_coordinate = (dim >= 2 ? 1 : 0);
2567  constexpr unsigned int z_coordinate = (dim >= 3 ? 2 : 0);
2568 
2569  switch (this->kind)
2570  {
2572  {
2573  // Vertices are special cases in that they do not actually
2574  // have coordinates. Error out if this function is called
2575  // with a vertex:
2576  Assert(false,
2577  ExcMessage("Vertices are zero-dimensional objects and "
2578  "as a consequence have no coordinates. You "
2579  "cannot meaningfully ask whether a point is "
2580  "inside a vertex (within a certain tolerance) "
2581  "without coordinate values."));
2582  return false;
2583  }
2584  case ReferenceCells::Line:
2587  {
2588  for (unsigned int d = 0; d < dim; ++d)
2589  if ((p[d] < -tolerance) || (p[d] > 1 + tolerance))
2590  return false;
2591  return true;
2592  }
2595  {
2596  // First make sure that we are in the first quadrant or octant
2597  for (unsigned int d = 0; d < dim; ++d)
2598  if (p[d] < -tolerance)
2599  return false;
2600 
2601  // Now we also need to make sure that we are below the diagonal line
2602  // or plane that delineates the simplex. This diagonal is given by
2603  // sum(p[d])<=1, and a diagonal a distance eps away is given by
2604  // sum(p[d])<=1+eps*sqrt(d). (For example, the point at (1,1) is a
2605  // distance of 1/sqrt(2) away from the diagonal. That is, its
2606  // sum satisfies
2607  // sum(p[d]) = 2 <= 1 + (1/sqrt(2)) * sqrt(2)
2608  // in other words, it satisfies the predicate with eps=1/sqrt(2).)
2609  double sum = 0;
2610  for (unsigned int d = 0; d < dim; ++d)
2611  sum += p[d];
2612  return (sum <= 1 + tolerance * std::sqrt(1. * dim));
2613  }
2615  {
2616  // A pyramid only lives in the upper half-space:
2617  if (p[z_coordinate] < -tolerance)
2618  return false;
2619 
2620  // It also only lives in the space below z=1:
2621  if (p[z_coordinate] > 1 + tolerance)
2622  return false;
2623 
2624  // Within what's left of the space, a pyramid is a cone that tapers
2625  // towards the top. First compute the distance of the point to the
2626  // axis in the max norm (this is the right norm because the vertices
2627  // of the pyramid are at points +/-1, +/-1):
2628  const double distance_from_axis =
2629  std::max(std::fabs(p[x_coordinate]), std::fabs(p[y_coordinate]));
2630 
2631  // We are inside the pyramid if the distance from the axis is less
2632  // than (1-z)
2633  return (distance_from_axis <= 1 + tolerance - p[z_coordinate]);
2634  }
2635  case ReferenceCells::Wedge:
2636  {
2637  // The wedge we use is a triangle extruded into the third
2638  // dimension by one unit. So we can use the same logic as for
2639  // triangles above (i.e., for the simplex above, using dim==2)
2640  // and then check the third dimension separately.
2641 
2642  if ((p[x_coordinate] < -tolerance) || (p[y_coordinate] < -tolerance))
2643  return false;
2644 
2645  const double sum = p[x_coordinate] + p[y_coordinate];
2646  if (sum > 1 + tolerance * std::sqrt(2.0))
2647  return false;
2648 
2649  if (p[z_coordinate] < -tolerance)
2650  return false;
2651  if (p[z_coordinate] > 1 + tolerance)
2652  return false;
2653 
2654  return true;
2655  }
2656  default:
2657  Assert(false, ExcNotImplemented());
2658  }
2659 
2660  return false;
2661 }
2662 
2663 
2664 
2665 template <int dim>
2666 inline Tensor<1, dim>
2667 ReferenceCell::unit_tangential_vectors(const unsigned int face_no,
2668  const unsigned int i) const
2669 {
2671  AssertIndexRange(i, dim - 1);
2672 
2673  switch (this->kind)
2674  {
2676  case ReferenceCells::Line:
2680  return GeometryInfo<dim>::unit_tangential_vectors[face_no][i];
2682  {
2683  AssertIndexRange(face_no, 3);
2684  static const std::array<Tensor<1, dim>, 3> table = {
2685  {Point<dim>(1, 0),
2686  Point<dim>(-std::sqrt(0.5), +std::sqrt(0.5)),
2687  Point<dim>(0, -1)}};
2688 
2689  return table[face_no];
2690  }
2692  {
2693  AssertIndexRange(face_no, 4);
2694  static const ndarray<Tensor<1, dim>, 4, 2> table = {
2695  {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
2696  {{Point<dim>(1, 0, 0), Point<dim>(0, 0, 1)}},
2697  {{Point<dim>(0, 0, 1), Point<dim>(0, 1, 0)}},
2698  {{Point<dim>(-std::pow(1.0 / 3.0, 1.0 / 4.0),
2699  +std::pow(1.0 / 3.0, 1.0 / 4.0),
2700  0),
2701  Point<dim>(-std::pow(1.0 / 3.0, 1.0 / 4.0),
2702  0,
2703  +std::pow(1.0 / 3.0, 1.0 / 4.0))}}}};
2704 
2705  return table[face_no][i];
2706  }
2708  {
2709  AssertIndexRange(face_no, 5);
2710  static const ndarray<Tensor<1, dim>, 5, 2> table = {
2711  {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
2712  {{Point<dim>(+1.0 / sqrt(2.0), 0, +1.0 / sqrt(2.0)),
2713  Point<dim>(0, 1, 0)}},
2714  {{Point<dim>(+1.0 / sqrt(2.0), 0, -1.0 / sqrt(2.0)),
2715  Point<dim>(0, 1, 0)}},
2716  {{Point<dim>(1, 0, 0),
2717  Point<dim>(0, +1.0 / sqrt(2.0), +1.0 / sqrt(2.0))}},
2718  {{Point<dim>(1, 0, 0),
2719  Point<dim>(0, +1.0 / sqrt(2.0), -1.0 / sqrt(2.0))}}}};
2720 
2721  return table[face_no][i];
2722  }
2723  case ReferenceCells::Wedge:
2724  {
2725  AssertIndexRange(face_no, 5);
2726  static const ndarray<Tensor<1, dim>, 5, 2> table = {
2727  {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
2728  {{Point<dim>(1, 0, 0), Point<dim>(0, 1, 0)}},
2729  {{Point<dim>(1, 0, 0), Point<dim>(0, 0, 1)}},
2730  {{Point<dim>(-1 / std::sqrt(2.0), +1 / std::sqrt(2.0), 0),
2731  Point<dim>(0, 0, 1)}},
2732  {{Point<dim>(0, 0, 1), Point<dim>(0, 1, 0)}}}};
2733 
2734  return table[face_no][i];
2735  }
2736  default:
2737  Assert(false, ExcNotImplemented());
2738  }
2739 
2740 
2741  return {};
2742 }
2743 
2744 
2745 
2746 template <int dim>
2747 inline Tensor<1, dim>
2748 ReferenceCell::unit_normal_vectors(const unsigned int face_no) const
2749 {
2750  AssertDimension(dim, this->get_dimension());
2751 
2752  if (is_hyper_cube())
2753  {
2755  return GeometryInfo<dim>::unit_normal_vector[face_no];
2756  }
2757  else if (dim == 2)
2758  {
2760 
2761  // Return the rotated vector
2762  return cross_product_2d(unit_tangential_vectors<dim>(face_no, 0));
2763  }
2764  else if (dim == 3)
2765  {
2766  return cross_product_3d(unit_tangential_vectors<dim>(face_no, 0),
2767  unit_tangential_vectors<dim>(face_no, 1));
2768  }
2769 
2770  Assert(false, ExcNotImplemented());
2771 
2772  return {};
2773 }
2774 
2775 
2776 
2777 inline unsigned int
2778 ReferenceCell::n_face_orientations(const unsigned int face_no) const
2779 {
2780  AssertIndexRange(face_no, n_faces());
2781  if (get_dimension() == 1)
2782  return 1;
2783  if (get_dimension() == 2)
2784  return 2;
2786  return 8;
2787  else if (face_reference_cell(face_no) == ReferenceCells::Triangle)
2788  return 6;
2789 
2790  Assert(false, ExcInternalError());
2792 }
2793 
2794 
2795 
2796 inline bool
2798  const unsigned int line,
2799  const unsigned char combined_face_orientation,
2800  const bool line_orientation) const
2801 {
2802  if (*this == ReferenceCells::Hexahedron)
2803  {
2804  static constexpr ::ndarray<bool, 2, 8> bool_table{
2805  {{{true, true, false, true, false, false, true, false}},
2806  {{true, true, true, false, false, false, false, true}}}};
2807 
2808  return (line_orientation ==
2809  bool_table[line / 2][combined_face_orientation]);
2810  }
2811  else
2812  // TODO: This might actually be wrong for some of the other
2813  // kinds of objects. We should check this
2814  return true;
2815 }
2816 
2817 
2818 
2819 namespace internal
2820 {
2821  template <typename T>
2823  {
2824  public:
2834  {
2837  }
2838 
2842  virtual ~NoPermutation() noexcept override = default;
2843 
2847  virtual void
2848  print_info(std::ostream &out) const override
2849  {
2850  out << '[';
2851 
2852  const unsigned int n_vertices = entity_type.n_vertices();
2853 
2854  for (unsigned int i = 0; i < n_vertices; ++i)
2855  {
2856  out << vertices_0[i];
2857  if (i + 1 != n_vertices)
2858  out << ',';
2859  }
2860 
2861  out << "] is not a valid permutation of [";
2862 
2863  for (unsigned int i = 0; i < n_vertices; ++i)
2864  {
2865  out << vertices_1[i];
2866  if (i + 1 != n_vertices)
2867  out << ',';
2868  }
2869 
2870  out << "]." << std::endl;
2871  }
2872 
2877 
2882 
2887  };
2888 } // namespace internal
2889 
2890 
2891 
2892 template <typename T, std::size_t N>
2893 inline unsigned char
2894 ReferenceCell::compute_orientation(const std::array<T, N> &vertices_0,
2895  const std::array<T, N> &vertices_1) const
2896 {
2897  Assert(N >= n_vertices(),
2898  ExcMessage("The number of array elements must be equal to or "
2899  "greater than the number of vertices of the cell "
2900  "referenced by this object."));
2901 
2902  // Call the non-deprecated function, taking care of calling it only with
2903  // those array elements that we actually care about (see the note
2904  // in the documentation about the arguments potentially being
2905  // larger arrays than necessary).
2906  return get_combined_orientation(
2907  make_array_view(vertices_0.begin(), vertices_0.begin() + n_vertices()),
2908  make_array_view(vertices_1.begin(), vertices_1.begin() + n_vertices()));
2909 }
2910 
2911 
2912 
2913 template <typename T>
2914 unsigned char
2916  const ArrayView<const T> &vertices_0,
2917  const ArrayView<const T> &vertices_1) const
2918 {
2919  Assert(vertices_0.size() == n_vertices(),
2920  ExcMessage("The number of array elements must be equal to "
2921  "the number of vertices of the cell "
2922  "referenced by this object."));
2923  Assert(vertices_1.size() == n_vertices(),
2924  ExcMessage("The number of array elements must be equal to "
2925  "the number of vertices of the cell "
2926  "referenced by this object."));
2927 
2928  const auto v0_equals = [&](const std::initializer_list<const T> &list) {
2929  Assert(list.size() == n_vertices(), ExcInternalError());
2930  return std::equal(vertices_0.begin(), vertices_0.end(), std::begin(list));
2931  };
2932 
2933  switch (this->kind)
2934  {
2935  case ReferenceCells::Line:
2936  // line_orientation=true
2937  if (v0_equals({vertices_1[0], vertices_1[1]}))
2938  return 1;
2939 
2940  // line_orientation=false
2941  if (v0_equals({vertices_1[1], vertices_1[0]}))
2942  return 0;
2943  break;
2945  // face_orientation=true, face_rotation=false, face_flip=false
2946  if (v0_equals({vertices_1[0], vertices_1[1], vertices_1[2]}))
2947  return 1;
2948 
2949  // face_orientation=true, face_rotation=true, face_flip=false
2950  if (v0_equals({vertices_1[1], vertices_1[2], vertices_1[0]}))
2951  return 3;
2952 
2953  // face_orientation=true, face_rotation=false, face_flip=true
2954  if (v0_equals({vertices_1[2], vertices_1[0], vertices_1[1]}))
2955  return 5;
2956 
2957  // face_orientation=false, face_rotation=false, face_flip=false
2958  if (v0_equals({vertices_1[0], vertices_1[2], vertices_1[1]}))
2959  return 0;
2960 
2961  // face_orientation=false, face_rotation=true, face_flip=false
2962  if (v0_equals({vertices_1[2], vertices_1[1], vertices_1[0]}))
2963  return 2;
2964 
2965  // face_orientation=false, face_rotation=false, face_flip=true
2966  if (v0_equals({vertices_1[1], vertices_1[0], vertices_1[2]}))
2967  return 4;
2968  break;
2970  // face_orientation=true, face_rotation=false, face_flip=false
2971  if (v0_equals(
2972  {vertices_1[0], vertices_1[1], vertices_1[2], vertices_1[3]}))
2973  return 1;
2974 
2975  // face_orientation=true, face_rotation=true, face_flip=false
2976  if (v0_equals(
2977  {vertices_1[2], vertices_1[0], vertices_1[3], vertices_1[1]}))
2978  return 3;
2979 
2980  // face_orientation=true, face_rotation=false, face_flip=true
2981  if (v0_equals(
2982  {vertices_1[3], vertices_1[2], vertices_1[1], vertices_1[0]}))
2983  return 5;
2984 
2985  // face_orientation=true, face_rotation=true, face_flip=true
2986  if (v0_equals(
2987  {vertices_1[1], vertices_1[3], vertices_1[0], vertices_1[2]}))
2988  return 7;
2989 
2990  // face_orientation=false, face_rotation=false, face_flip=false
2991  if (v0_equals(
2992  {vertices_1[0], vertices_1[2], vertices_1[1], vertices_1[3]}))
2993  return 0;
2994 
2995  // face_orientation=false, face_rotation=true, face_flip=false
2996  if (v0_equals(
2997  {vertices_1[2], vertices_1[3], vertices_1[0], vertices_1[1]}))
2998  return 2;
2999 
3000  // face_orientation=false, face_rotation=false, face_flip=true
3001  if (v0_equals(
3002  {vertices_1[3], vertices_1[1], vertices_1[2], vertices_1[0]}))
3003  return 4;
3004 
3005  // face_orientation=false, face_rotation=true, face_flip=true
3006  if (v0_equals(
3007  {vertices_1[1], vertices_1[0], vertices_1[3], vertices_1[2]}))
3008  return 6;
3009  break;
3010  default:
3011  Assert(false, ExcNotImplemented());
3012  }
3013 
3014  Assert(false, (internal::NoPermutation<T>(*this, vertices_0, vertices_1)));
3016 }
3017 
3018 
3019 
3020 template <typename T, std::size_t N>
3021 inline std::array<T, N>
3023  const std::array<T, N> &vertices,
3024  const unsigned int orientation) const
3025 {
3026  Assert(N >= n_vertices(),
3027  ExcMessage("The number of array elements must be equal to or "
3028  "greater than the number of vertices of the cell "
3029  "referenced by this object."));
3030 
3031  // Call the non-deprecated function, taking care of calling it only with
3032  // those array elements that we actually care about (see the note
3033  // in the documentation about the arguments potentially being
3034  // larger arrays than necessary).
3035  const auto permutation = permute_by_combined_orientation(
3036  make_array_view(vertices.begin(), vertices.begin() + n_vertices()),
3037  orientation);
3038 
3039  std::array<T, N> temp;
3040  std::copy(permutation.begin(), permutation.end(), temp.begin());
3041 
3042  return temp;
3043 }
3044 
3045 
3046 
3047 template <typename T>
3048 boost::container::small_vector<T, 8>
3051  const unsigned char orientation) const
3052 {
3053  Assert(vertices.size() == n_vertices(),
3054  ExcMessage("The number of array elements must be equal to "
3055  "the number of vertices of the cell "
3056  "referenced by this object."));
3057 
3058  switch (this->kind)
3059  {
3060  case ReferenceCells::Line:
3061  switch (orientation)
3062  {
3063  case 1:
3064  return {vertices[0], vertices[1]};
3065  case 0:
3066  return {vertices[1], vertices[0]};
3067  default:
3068  Assert(false, ExcNotImplemented());
3069  }
3070  break;
3072  switch (orientation)
3073  {
3074  case 1:
3075  return {vertices[0], vertices[1], vertices[2]};
3076  case 3:
3077  return {vertices[1], vertices[2], vertices[0]};
3078  case 5:
3079  return {vertices[2], vertices[0], vertices[1]};
3080  case 0:
3081  return {vertices[0], vertices[2], vertices[1]};
3082  case 2:
3083  return {vertices[2], vertices[1], vertices[0]};
3084  case 4:
3085  return {vertices[1], vertices[0], vertices[2]};
3086  default:
3087  Assert(false, ExcNotImplemented());
3088  }
3089  break;
3091  switch (orientation)
3092  {
3093  case 1:
3094  return {vertices[0], vertices[1], vertices[2], vertices[3]};
3095  case 3:
3096  return {vertices[2], vertices[0], vertices[3], vertices[1]};
3097  case 5:
3098  return {vertices[3], vertices[2], vertices[1], vertices[0]};
3099  case 7:
3100  return {vertices[1], vertices[3], vertices[0], vertices[2]};
3101  case 0:
3102  return {vertices[0], vertices[2], vertices[1], vertices[3]};
3103  case 2:
3104  return {vertices[2], vertices[3], vertices[0], vertices[1]};
3105  case 4:
3106  return {vertices[3], vertices[1], vertices[2], vertices[0]};
3107  case 6:
3108  return {vertices[1], vertices[0], vertices[3], vertices[2]};
3109  default:
3110  Assert(false, ExcNotImplemented());
3111  }
3112  break;
3113  default:
3114  AssertThrow(false, ExcNotImplemented());
3115  }
3116 
3117  return {};
3118 }
3119 
3120 
3121 
3122 template <>
3123 unsigned int
3124 ReferenceCell::vtk_lexicographic_to_node_index<0>(
3125  const std::array<unsigned, 0> &node_indices,
3126  const std::array<unsigned, 0> &nodes_per_direction,
3127  const bool legacy_format) const;
3128 
3129 template <>
3130 unsigned int
3131 ReferenceCell::vtk_lexicographic_to_node_index<1>(
3132  const std::array<unsigned, 1> &node_indices,
3133  const std::array<unsigned, 1> &nodes_per_direction,
3134  const bool legacy_format) const;
3135 
3136 template <>
3137 unsigned int
3138 ReferenceCell::vtk_lexicographic_to_node_index<2>(
3139  const std::array<unsigned, 2> &node_indices,
3140  const std::array<unsigned, 2> &nodes_per_direction,
3141  const bool legacy_format) const;
3142 
3143 template <>
3144 unsigned int
3145 ReferenceCell::vtk_lexicographic_to_node_index<3>(
3146  const std::array<unsigned, 3> &node_indices,
3147  const std::array<unsigned, 3> &nodes_per_direction,
3148  const bool legacy_format) const;
3149 
3151 
3152 #endif
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Definition: array_view.h:704
OutputOperator< VectorType > & operator<<(OutputOperator< VectorType > &out, unsigned int step)
Definition: operator.h:165
iterator begin() const
Definition: array_view.h:594
iterator end() const
Definition: array_view.h:603
std::size_t size() const
Definition: array_view.h:576
Abstract base class for mapping classes.
Definition: mapping.h:317
Definition: point.h:112
std::istream & operator>>(std::istream &in, Point< dim, Number > &p)
Definition: point.h:699
static Point< dim, Number > unit_vector(const unsigned int i)
ArrayView< const unsigned int > faces_for_given_vertex(const unsigned int vertex_index) const
unsigned int child_cell_on_face(const unsigned int face, const unsigned int subface, const unsigned char face_orientation=default_combined_face_orientation()) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
unsigned int n_vertices() const
boost::container::small_vector< T, 8 > permute_by_combined_orientation(const ArrayView< const T > &vertices, const unsigned char orientation) const
Quadrature< dim > get_gauss_type_quadrature(const unsigned n_points_1d) const
bool is_hyper_cube() const
unsigned int vtk_linear_type() const
friend std::ostream & operator<<(std::ostream &out, const ReferenceCell &reference_cell)
void serialize(Archive &archive, const unsigned int)
friend std::istream & operator>>(std::istream &in, ReferenceCell &reference_cell)
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const
static constexpr unsigned char default_combined_face_orientation()
unsigned char compute_orientation(const std::array< T, N > &vertices_0, const std::array< T, N > &vertices_1) const
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const
bool standard_vs_true_line_orientation(const unsigned int line, const unsigned char face_orientation, const bool line_orientation) const
Point< dim > vertex(const unsigned int v) const
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char face_orientation) const
const Quadrature< dim > & get_nodal_type_quadrature() const
Point< dim > face_vertex_location(const unsigned int face, const unsigned int vertex) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > isotropic_child_indices() const
std::array< unsigned int, 2 > standard_line_to_face_and_line_index(const unsigned int line) const
double d_linear_shape_function(const Point< dim > &xi, const unsigned int i) const
Tensor< 1, dim > unit_tangential_vectors(const unsigned int face_no, const unsigned int i) const
static constexpr ndarray< unsigned int, 2, 2 > line_vertex_permutations
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
std::array< T, N > permute_according_orientation(const std::array< T, N > &vertices, const unsigned int orientation) const
bool contains_point(const Point< dim > &p, const double tolerance=0) const
unsigned int n_isotropic_children() const
unsigned int gmsh_element_type() const
unsigned int n_face_orientations(const unsigned int face_no) const
static constexpr ndarray< unsigned int, 8, 4 > quadrilateral_vertex_permutations
double volume() const
unsigned int n_faces() const
unsigned int line_to_cell_vertices(const unsigned int line, const unsigned int vertex) const
std::uint8_t kind
constexpr ReferenceCell()
unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const
unsigned int unv_vertex_to_deal_vertex(const unsigned int vertex_n) const
unsigned int vtk_vertex_to_deal_vertex(const unsigned int vertex_index) const
std::unique_ptr< Mapping< dim, spacedim > > get_default_mapping(const unsigned int degree) const
Quadrature< dim > get_midpoint_quadrature() const
unsigned int vtk_quadratic_type() const
static constexpr unsigned char reversed_combined_line_orientation()
unsigned int n_lines() const
unsigned int vtk_lagrange_type() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
unsigned int get_dimension() const
ReferenceCell face_reference_cell(const unsigned int face_no) const
unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const
static ReferenceCell n_vertices_to_type(const int dim, const unsigned int n_vertices)
unsigned int standard_to_real_face_line(const unsigned int line, const unsigned int face, const unsigned char face_orientation) const
static constexpr std::size_t memory_consumption()
const Mapping< dim, spacedim > & get_default_linear_mapping() const
std::string to_string() const
bool is_simplex() const
unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const unsigned char face_orientation) const
constexpr bool operator!=(const ReferenceCell &type) const
unsigned char get_combined_orientation(const ArrayView< const T > &vertices_0, const ArrayView< const T > &vertices_1) const
constexpr bool operator==(const ReferenceCell &type) const
Point< dim > barycenter() const
Tensor< 1, dim > unit_normal_vectors(const unsigned int face_no) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
Point< dim > closest_point(const Point< dim > &p) const
static constexpr ndarray< unsigned int, 6, 3 > triangle_vertex_permutations
Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i) const
virtual void print_info(std::ostream &out) const override
const ReferenceCell entity_type
virtual ~NoPermutation() noexcept override=default
const ArrayView< const T > vertices_0
const ArrayView< const T > vertices_1
NoPermutation(const ReferenceCell &entity_type, const ArrayView< const T > &vertices_0, const ArrayView< const T > &vertices_1)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:474
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
Point< 3 > vertices[4]
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1614
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1787
#define AssertIndexRange(index, range)
Definition: exceptions.h:1855
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1703
Expression fabs(const Expression &x)
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
static const char U
static const char N
std::string to_string(const T &t)
Definition: patterns.h:2392
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
constexpr const ReferenceCell Tetrahedron
constexpr const ReferenceCell Quadrilateral
constexpr const ReferenceCell Wedge
constexpr const ReferenceCell Pyramid
constexpr const ReferenceCell Invalid
constexpr const ReferenceCell Triangle
constexpr const ReferenceCell Hexahedron
constexpr const ReferenceCell & get_hypercube()
constexpr const ReferenceCell Vertex
constexpr const ReferenceCell Line
constexpr const ReferenceCell & get_simplex()
VectorType::value_type * begin(VectorType &V)
T sum(const T &t, const MPI_Comm mpi_communicator)
bool get_bit(const unsigned char number, const unsigned int n)
Definition: utilities.h:1566
void copy(const T *begin, const T *end, U *dest)
constexpr ReferenceCell make_reference_cell_from_int(const std::uint8_t kind)
static const unsigned int invalid_unsigned_int
Definition: types.h:213
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
typename internal::ndarray::HelperArray< T, Ns... >::type ndarray
Definition: ndarray.h:108
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)
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)
static std::array< unsigned int, 2 > standard_hex_line_to_quad_line_index(const unsigned int line)
static unsigned int standard_to_real_face_line(const unsigned int line, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
static double d_linear_shape_function(const Point< dim > &xi, const unsigned int i)
static Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i)
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)
static std::array< unsigned int, 2 > standard_quad_vertex_to_line_vertex_index(const unsigned int vertex)
static std::array< unsigned int, 2 > standard_hex_vertex_to_quad_vertex_index(const unsigned int vertex)