Reference documentation for deal.II version Git 62a7c2142d 2021-03-08 22:48:39 +0100
\(\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 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/tensor.h>
25 #include <deal.II/base/utilities.h>
26 
27 
29 
30 // Forward declarations
31 #ifndef DOXYGEN
32 template <int dim, int spacedim>
33 class Mapping;
34 
35 template <int dim>
36 class Quadrature;
37 
38 class ReferenceCell;
39 #endif
40 
41 
42 namespace internal
43 {
44  namespace ReferenceCell
45  {
59  DEAL_II_CONSTEXPR ::ReferenceCell
60  make_reference_cell_from_int(const std::uint8_t kind);
61  } // namespace ReferenceCell
62 } // namespace internal
63 
64 
65 
87 {
88 public:
96  static ReferenceCell
97  n_vertices_to_type(const int dim, const unsigned int n_vertices);
98 
108  ReferenceCell();
109 
120  bool
121  is_hyper_cube() const;
122 
126  bool
127  is_simplex() const;
128 
133  unsigned int
134  get_dimension() const;
135 
149  template <int dim>
150  double
151  d_linear_shape_function(const Point<dim> &xi, const unsigned int i) const;
152 
157  template <int dim>
159  d_linear_shape_function_gradient(const Point<dim> & xi,
160  const unsigned int i) const;
161 
171  template <int dim, int spacedim>
172  std::unique_ptr<Mapping<dim, spacedim>>
173  get_default_mapping(const unsigned int degree) const;
174 
186  template <int dim, int spacedim>
187  const Mapping<dim, spacedim> &
189 
198  template <int dim>
200  get_gauss_type_quadrature(const unsigned n_points_1D) const;
201 
208  template <int dim>
209  const Quadrature<dim> &
210  get_nodal_type_quadrature() const;
211 
226  unsigned int
227  n_vertices() const;
228 
234  vertex_indices() const;
235 
241  unsigned int
242  n_lines() const;
243 
249  line_indices() const;
250 
256  unsigned int
257  n_faces() const;
258 
264  face_indices() const;
265 
278  ReferenceCell
279  face_reference_cell(const unsigned int face_no) const;
280 
322  unsigned int
323  child_cell_on_face(const unsigned int face_n,
324  const unsigned int subface_n) const;
325 
335  std::array<unsigned int, 2>
336  standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const;
337 
347  std::array<unsigned int, 2>
348  standard_line_to_face_and_line_index(const unsigned int line) const;
349 
353  unsigned int
354  face_to_cell_lines(const unsigned int face,
355  const unsigned int line,
356  const unsigned char face_orientation) const;
357 
361  unsigned int
362  face_to_cell_vertices(const unsigned int face,
363  const unsigned int vertex,
364  const unsigned char face_orientation) const;
365 
369  unsigned int
370  standard_to_real_face_vertex(const unsigned int vertex,
371  const unsigned int face,
372  const unsigned char face_orientation) const;
373 
377  unsigned int
378  standard_to_real_face_line(const unsigned int line,
379  const unsigned int face,
380  const unsigned char face_orientation) const;
381 
392  bool
393  standard_vs_true_line_orientation(const unsigned int line,
394  const unsigned char face_orientation,
395  const unsigned char line_orientation) const;
396 
407  /*
408  * Return @f$i@f$-th unit tangential vector of a face of the reference cell.
409  * The vectors are arranged such that the
410  * cross product between the two vectors returns the unit normal vector.
411  *
412  * @pre @f$i@f$ must be between zero and `dim-1`.
413  */
414  template <int dim>
416  unit_tangential_vectors(const unsigned int face_no,
417  const unsigned int i) const;
418 
422  template <int dim>
424  unit_normal_vectors(const unsigned int face_no) const;
425 
430  template <typename T, std::size_t N>
431  unsigned char
432  compute_orientation(const std::array<T, N> &vertices_0,
433  const std::array<T, N> &vertices_1) const;
434 
438  template <typename T, std::size_t N>
439  std::array<T, N>
440  permute_according_orientation(const std::array<T, N> &vertices,
441  const unsigned int orientation) const;
442 
447  faces_for_given_vertex(const unsigned int vertex_index) const;
448 
461  unsigned int
462  exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const;
463 
467  unsigned int
468  exodusii_face_to_deal_face(const unsigned int face_n) const;
469 
483  std::string
484  to_string() const;
485 
489  constexpr operator std::uint8_t() const;
490 
494  constexpr bool
495  operator==(const ReferenceCell &type) const;
496 
500  constexpr bool
501  operator!=(const ReferenceCell &type) const;
502 
508  template <class Archive>
509  void
510  serialize(Archive &archive, const unsigned int /*version*/);
511 
516 private:
520  std::uint8_t kind;
521 
528  constexpr ReferenceCell(const std::uint8_t kind);
529 
534  friend DEAL_II_CONSTEXPR ReferenceCell
536 };
537 
538 
539 
540 inline constexpr ReferenceCell::ReferenceCell(const std::uint8_t kind)
541  : kind(kind)
542 {}
543 
544 
545 
546 inline constexpr ReferenceCell::operator std::uint8_t() const
547 {
548  return kind;
549 }
550 
551 
552 
553 inline constexpr bool
555 {
556  return kind == type.kind;
557 }
558 
559 
560 
561 inline constexpr bool
563 {
564  return kind != type.kind;
565 }
566 
567 
568 
569 namespace internal
570 {
571  namespace ReferenceCell
572  {
573  inline DEAL_II_CONSTEXPR ::ReferenceCell
574  make_reference_cell_from_int(const std::uint8_t kind)
575  {
576  // Make sure these are the only indices from which objects can be
577  // created.
578  Assert((kind == static_cast<std::uint8_t>(-1)) || (kind < 8),
579  ExcInternalError());
580 
581  // Call the private constructor, which we can from here because this
582  // function is a 'friend'.
583  return {kind};
584  }
585  } // namespace ReferenceCell
586 } // namespace internal
587 
588 
589 
597 namespace ReferenceCells
598 {
617  static_cast<std::uint8_t>(-1));
618 
624  template <int dim>
625  constexpr const ReferenceCell &
626  get_simplex();
627 
633  template <int dim>
634  constexpr const ReferenceCell &
635  get_hypercube();
636 } // namespace ReferenceCells
637 
638 
639 
640 inline DEAL_II_CONSTEXPR
642  : ReferenceCell(ReferenceCells::Invalid)
643 {}
644 
645 
646 
647 template <class Archive>
648 inline void
649 ReferenceCell::serialize(Archive &archive, const unsigned int /*version*/)
650 {
651  archive &kind;
652 }
653 
654 
655 
657 ReferenceCell::faces_for_given_vertex(const unsigned int vertex) const
658 {
659  if (*this == ReferenceCells::Line)
660  {
662  return {&GeometryInfo<2>::vertex_to_face[vertex][0], 1};
663  }
664  else if (*this == ReferenceCells::Quadrilateral)
665  {
667  return {&GeometryInfo<2>::vertex_to_face[vertex][0], 2};
668  }
669  else if (*this == ReferenceCells::Hexahedron)
670  {
672  return {&GeometryInfo<3>::vertex_to_face[vertex][0], 3};
673  }
674  else if (*this == ReferenceCells::Triangle)
675  {
676  AssertIndexRange(vertex, 3);
677  static const ndarray<unsigned int, 3, 2> table = {
678  {{{0, 2}}, {{0, 1}}, {{1, 2}}}};
679 
680  return table[vertex];
681  }
682  else if (*this == ReferenceCells::Tetrahedron)
683  {
684  AssertIndexRange(vertex, 4);
685  static const ndarray<unsigned int, 4, 3> table = {
686  {{{0, 1, 2}}, {{0, 1, 3}}, {{0, 2, 3}}, {{1, 2, 3}}}};
687 
688  return table[vertex];
689  }
690  else if (*this == ReferenceCells::Wedge)
691  {
692  AssertIndexRange(vertex, 6);
693  static const ndarray<unsigned int, 6, 3> table = {{{{0, 2, 4}},
694  {{0, 2, 3}},
695  {{0, 3, 4}},
696  {{1, 2, 4}},
697  {{1, 2, 3}},
698  {{1, 3, 4}}}};
699 
700  return table[vertex];
701  }
702  else if (*this == ReferenceCells::Pyramid)
703  {
704  AssertIndexRange(vertex, 5);
705  static const unsigned int X = numbers::invalid_unsigned_int;
706  static const ndarray<unsigned int, 5, 4> table = {{{{0, 1, 3, X}},
707  {{0, 2, 3, X}},
708  {{0, 1, 4, X}},
709  {{0, 2, 4, X}},
710  {{1, 2, 3, 4}}}};
711 
712  return {&table[vertex][0], vertex == 4 ? 4u : 3u};
713  }
714 
715  Assert(false, ExcNotImplemented());
716 
717  return {};
718 }
719 
720 
721 
722 inline bool
724 {
725  return (*this == ReferenceCells::Vertex || *this == ReferenceCells::Line ||
727  *this == ReferenceCells::Hexahedron);
728 }
729 
730 
731 
732 inline bool
734 {
735  return (*this == ReferenceCells::Vertex || *this == ReferenceCells::Line ||
736  *this == ReferenceCells::Triangle ||
737  *this == ReferenceCells::Tetrahedron);
738 }
739 
740 
741 
742 inline unsigned int
744 {
745  if (*this == ReferenceCells::Vertex)
746  return 0;
747  else if (*this == ReferenceCells::Line)
748  return 1;
749  else if ((*this == ReferenceCells::Triangle) ||
751  return 2;
752  else if ((*this == ReferenceCells::Tetrahedron) ||
753  (*this == ReferenceCells::Pyramid) ||
754  (*this == ReferenceCells::Wedge) ||
755  (*this == ReferenceCells::Hexahedron))
756  return 3;
757 
758  Assert(false, ExcNotImplemented());
760 }
761 
762 
763 
764 inline unsigned int
766 {
767  if (*this == ReferenceCells::Vertex)
768  return 1;
769  else if (*this == ReferenceCells::Line)
770  return 2;
771  else if (*this == ReferenceCells::Triangle)
772  return 3;
773  else if (*this == ReferenceCells::Quadrilateral)
774  return 4;
775  else if (*this == ReferenceCells::Tetrahedron)
776  return 4;
777  else if (*this == ReferenceCells::Pyramid)
778  return 5;
779  else if (*this == ReferenceCells::Wedge)
780  return 6;
781  else if (*this == ReferenceCells::Hexahedron)
782  return 8;
783 
784  Assert(false, ExcNotImplemented());
786 }
787 
788 
789 
790 inline unsigned int
792 {
793  if (*this == ReferenceCells::Vertex)
794  return 0;
795  else if (*this == ReferenceCells::Line)
796  return 1;
797  else if (*this == ReferenceCells::Triangle)
798  return 3;
799  else if (*this == ReferenceCells::Quadrilateral)
800  return 4;
801  else if (*this == ReferenceCells::Tetrahedron)
802  return 6;
803  else if (*this == ReferenceCells::Pyramid)
804  return 7;
805  else if (*this == ReferenceCells::Wedge)
806  return 9;
807  else if (*this == ReferenceCells::Hexahedron)
808  return 12;
809 
810  Assert(false, ExcNotImplemented());
812 }
813 
814 
815 
816 inline unsigned int
818 {
819  if (*this == ReferenceCells::Vertex)
820  return 0;
821  else if (*this == ReferenceCells::Line)
822  return 2;
823  else if (*this == ReferenceCells::Triangle)
824  return 3;
825  else if (*this == ReferenceCells::Quadrilateral)
826  return 4;
827  else if (*this == ReferenceCells::Tetrahedron)
828  return 4;
829  else if (*this == ReferenceCells::Pyramid)
830  return 5;
831  else if (*this == ReferenceCells::Wedge)
832  return 5;
833  else if (*this == ReferenceCells::Hexahedron)
834  return 6;
835 
836  Assert(false, ExcNotImplemented());
838 }
839 
840 
841 
844 {
845  return {0U, n_vertices()};
846 }
847 
848 
849 
852 {
853  return {0U, n_lines()};
854 }
855 
856 
857 
860 {
861  return {0U, n_faces()};
862 }
863 
864 
865 
866 inline ReferenceCell
867 ReferenceCell::face_reference_cell(const unsigned int face_no) const
868 {
869  AssertIndexRange(face_no, n_faces());
870 
871  if (*this == ReferenceCells::Vertex)
873  else if (*this == ReferenceCells::Line)
874  return ReferenceCells::Vertex;
875  else if (*this == ReferenceCells::Triangle)
876  return ReferenceCells::Line;
877  else if (*this == ReferenceCells::Quadrilateral)
878  return ReferenceCells::Line;
879  else if (*this == ReferenceCells::Tetrahedron)
881  else if (*this == ReferenceCells::Pyramid)
882  {
883  if (face_no == 0)
885  else
887  }
888  else if (*this == ReferenceCells::Wedge)
889  {
890  if (face_no > 1)
892  else
894  }
895  else if (*this == ReferenceCells::Hexahedron)
897 
898  Assert(false, ExcNotImplemented());
900 }
901 
902 
903 
904 inline unsigned int
905 ReferenceCell::child_cell_on_face(const unsigned int face,
906  const unsigned int subface) const
907 {
908  AssertIndexRange(face, n_faces());
909 
910  if (*this == ReferenceCells::Vertex)
911  {
912  Assert(false, ExcNotImplemented());
913  }
914  else if (*this == ReferenceCells::Line)
915  {
916  Assert(false, ExcNotImplemented());
917  }
918  else if (*this == ReferenceCells::Triangle)
919  {
920  static const ndarray<unsigned int, 3, 2> subcells = {
921  {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
922 
923  return subcells[face][subface];
924  }
925  else if (*this == ReferenceCells::Quadrilateral)
926  {
927  Assert(false, ExcNotImplemented());
928  }
929  else if (*this == ReferenceCells::Tetrahedron)
930  {
931  Assert(false, ExcNotImplemented());
932  }
933  else if (*this == ReferenceCells::Pyramid)
934  {
935  Assert(false, ExcNotImplemented());
936  }
937  else if (*this == ReferenceCells::Wedge)
938  {
939  Assert(false, ExcNotImplemented());
940  }
941  else if (*this == ReferenceCells::Hexahedron)
942  {
943  Assert(false, ExcNotImplemented());
944  }
945 
946  Assert(false, ExcNotImplemented());
947  return {};
948 }
949 
950 
951 
952 inline std::array<unsigned int, 2>
954  const unsigned int vertex) const
955 {
956  AssertIndexRange(vertex, n_vertices());
957 
958  if (*this == ReferenceCells::Vertex)
959  {
960  Assert(false, ExcNotImplemented());
961  }
962  else if (*this == ReferenceCells::Line)
963  {
964  Assert(false, ExcNotImplemented());
965  }
966  else if (*this == ReferenceCells::Triangle)
967  {
968  static const ndarray<unsigned int, 3, 2> table = {
969  {{{0, 0}}, {{0, 1}}, {{1, 1}}}};
970 
971  return table[vertex];
972  }
973  else if (*this == ReferenceCells::Quadrilateral)
974  {
976  }
977  else if (*this == ReferenceCells::Tetrahedron)
978  {
979  static const ndarray<unsigned int, 4, 2> table = {
980  {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 2}}}};
981 
982  return table[vertex];
983  }
984  else if (*this == ReferenceCells::Pyramid)
985  {
986  static const ndarray<unsigned int, 5, 2> table = {
987  {{{0, 0}}, {{0, 1}}, {{0, 2}}, {{0, 3}}, {{1, 2}}}};
988 
989  return table[vertex];
990  }
991  else if (*this == ReferenceCells::Wedge)
992  {
993  static const ndarray<unsigned int, 6, 2> table = {
994  {{{0, 1}}, {{0, 0}}, {{0, 2}}, {{1, 0}}, {{1, 1}}, {{1, 2}}}};
995 
996  return table[vertex];
997  }
998  else if (*this == ReferenceCells::Hexahedron)
999  {
1001  }
1002 
1003  Assert(false, ExcNotImplemented());
1004  return {};
1005 }
1006 
1007 
1008 
1009 inline std::array<unsigned int, 2>
1011  const unsigned int line) const
1012 {
1013  AssertIndexRange(line, n_lines());
1014 
1015  if (*this == ReferenceCells::Vertex)
1016  {
1017  Assert(false, ExcNotImplemented());
1018  }
1019  else if (*this == ReferenceCells::Line)
1020  {
1021  Assert(false, ExcNotImplemented());
1022  }
1023  else if (*this == ReferenceCells::Triangle)
1024  {
1025  Assert(false, ExcNotImplemented());
1026  }
1027  else if (*this == ReferenceCells::Quadrilateral)
1028  {
1029  Assert(false, ExcNotImplemented());
1030  }
1031  else if (*this == ReferenceCells::Tetrahedron)
1032  {
1033  static const std::array<unsigned int, 2> table[6] = {
1034  {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}};
1035 
1036  return table[line];
1037  }
1038  else if (*this == ReferenceCells::Pyramid)
1039  {
1040  static const std::array<unsigned int, 2> table[8] = {{{0, 0}},
1041  {{0, 1}},
1042  {{0, 2}},
1043  {{0, 3}},
1044  {{1, 2}},
1045  {{2, 1}},
1046  {{1, 1}},
1047  {{2, 2}}};
1048 
1049  return table[line];
1050  }
1051  else if (*this == ReferenceCells::Wedge)
1052  {
1053  static const std::array<unsigned int, 2> table[9] = {{{0, 0}},
1054  {{0, 2}},
1055  {{0, 1}},
1056  {{1, 0}},
1057  {{1, 1}},
1058  {{1, 2}},
1059  {{2, 0}},
1060  {{2, 1}},
1061  {{3, 1}}};
1062 
1063  return table[line];
1064  }
1065  else if (*this == ReferenceCells::Hexahedron)
1066  {
1068  }
1069 
1070  Assert(false, ExcNotImplemented());
1071  return {};
1072 }
1073 
1074 
1075 
1076 inline unsigned int
1077 ReferenceCell::face_to_cell_lines(const unsigned int face,
1078  const unsigned int line,
1079  const unsigned char face_orientation) const
1080 {
1081  AssertIndexRange(face, n_faces());
1083 
1084  if (*this == ReferenceCells::Vertex)
1085  {
1086  Assert(false, ExcNotImplemented());
1087  }
1088  else if (*this == ReferenceCells::Line)
1089  {
1091  face,
1092  line,
1093  Utilities::get_bit(face_orientation, 0),
1094  Utilities::get_bit(face_orientation, 2),
1095  Utilities::get_bit(face_orientation, 1));
1096  }
1097  else if (*this == ReferenceCells::Triangle)
1098  {
1099  return face;
1100  }
1101  else if (*this == ReferenceCells::Quadrilateral)
1102  {
1104  face,
1105  line,
1106  Utilities::get_bit(face_orientation, 0),
1107  Utilities::get_bit(face_orientation, 2),
1108  Utilities::get_bit(face_orientation, 1));
1109  }
1110  else if (*this == ReferenceCells::Tetrahedron)
1111  {
1112  const static ndarray<unsigned int, 4, 3> table = {
1113  {{{0, 1, 2}}, {{0, 3, 4}}, {{2, 5, 3}}, {{1, 4, 5}}}};
1114 
1115  return table[face]
1116  [standard_to_real_face_line(line, face, face_orientation)];
1117  }
1118  else if (*this == ReferenceCells::Pyramid)
1119  {
1120  Assert(false, ExcNotImplemented());
1121  }
1122  else if (*this == ReferenceCells::Wedge)
1123  {
1124  Assert(false, ExcNotImplemented());
1125  }
1126  else if (*this == ReferenceCells::Hexahedron)
1127  {
1129  face,
1130  line,
1131  Utilities::get_bit(face_orientation, 0),
1132  Utilities::get_bit(face_orientation, 2),
1133  Utilities::get_bit(face_orientation, 1));
1134  }
1135 
1136  Assert(false, ExcNotImplemented());
1138 }
1139 
1140 
1141 
1142 inline unsigned int
1144  const unsigned int vertex,
1145  const unsigned char face_orientation) const
1146 {
1147  AssertIndexRange(face, n_faces());
1149 
1150  if (*this == ReferenceCells::Vertex)
1151  {
1152  Assert(false, ExcNotImplemented());
1153  }
1154  else if (*this == ReferenceCells::Line)
1155  {
1157  face,
1158  vertex,
1159  Utilities::get_bit(face_orientation, 0),
1160  Utilities::get_bit(face_orientation, 2),
1161  Utilities::get_bit(face_orientation, 1));
1162  }
1163  else if (*this == ReferenceCells::Triangle)
1164  {
1165  static const ndarray<unsigned int, 3, 2> table = {
1166  {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
1167 
1168  return table[face][face_orientation ? vertex : (1 - vertex)];
1169  }
1170  else if (*this == ReferenceCells::Quadrilateral)
1171  {
1173  face,
1174  vertex,
1175  Utilities::get_bit(face_orientation, 0),
1176  Utilities::get_bit(face_orientation, 2),
1177  Utilities::get_bit(face_orientation, 1));
1178  }
1179  else if (*this == ReferenceCells::Tetrahedron)
1180  {
1181  static const ndarray<unsigned int, 4, 3> table = {
1182  {{{0, 1, 2}}, {{1, 0, 3}}, {{0, 2, 3}}, {{2, 1, 3}}}};
1183 
1184  return table[face][standard_to_real_face_vertex(
1185  vertex, face, face_orientation)];
1186  }
1187  else if (*this == ReferenceCells::Pyramid)
1188  {
1189  constexpr auto X = numbers::invalid_unsigned_int;
1190  static const ndarray<unsigned int, 5, 4> table = {{{{0, 1, 2, 3}},
1191  {{0, 2, 4, X}},
1192  {{3, 1, 4, X}},
1193  {{1, 0, 4, X}},
1194  {{2, 3, 4, X}}}};
1195 
1196  return table[face][standard_to_real_face_vertex(
1197  vertex, face, face_orientation)];
1198  }
1199  else if (*this == ReferenceCells::Wedge)
1200  {
1201  constexpr auto X = numbers::invalid_unsigned_int;
1202  static const ndarray<unsigned int, 6, 4> table = {{{{1, 0, 2, X}},
1203  {{3, 4, 5, X}},
1204  {{0, 1, 3, 4}},
1205  {{1, 2, 4, 5}},
1206  {{2, 0, 5, 3}}}};
1207 
1208  return table[face][standard_to_real_face_vertex(
1209  vertex, face, face_orientation)];
1210  }
1211  else if (*this == ReferenceCells::Hexahedron)
1212  {
1214  face,
1215  vertex,
1216  Utilities::get_bit(face_orientation, 0),
1217  Utilities::get_bit(face_orientation, 2),
1218  Utilities::get_bit(face_orientation, 1));
1219  }
1220 
1221  Assert(false, ExcNotImplemented());
1223 }
1224 
1225 
1226 
1227 inline unsigned int
1229  const unsigned int vertex,
1230  const unsigned int face,
1231  const unsigned char face_orientation) const
1232 {
1233  AssertIndexRange(face, n_faces());
1235 
1236  if (*this == ReferenceCells::Vertex)
1237  {
1238  Assert(false, ExcNotImplemented());
1239  }
1240  else if (*this == ReferenceCells::Line)
1241  {
1242  Assert(false, ExcNotImplemented());
1243  }
1244  else if (*this == ReferenceCells::Triangle)
1245  {
1246  static const ndarray<unsigned int, 2, 2> table = {{{{1, 0}}, {{0, 1}}}};
1247 
1248  return table[face_orientation][vertex];
1249  }
1250  else if (*this == ReferenceCells::Quadrilateral)
1251  {
1253  face_orientation);
1254  }
1255  else if (*this == ReferenceCells::Tetrahedron)
1256  {
1257  static const ndarray<unsigned int, 6, 3> table = {{{{0, 2, 1}},
1258  {{0, 1, 2}},
1259  {{2, 1, 0}},
1260  {{1, 2, 0}},
1261  {{1, 0, 2}},
1262  {{2, 0, 1}}}};
1263 
1264  return table[face_orientation][vertex];
1265  }
1266  else if (*this == ReferenceCells::Pyramid)
1267  {
1268  if (face == 0) // The quadrilateral face
1269  {
1271  vertex,
1272  Utilities::get_bit(face_orientation, 0),
1273  Utilities::get_bit(face_orientation, 2),
1274  Utilities::get_bit(face_orientation, 1));
1275  }
1276  else // One of the triangular faces
1277  {
1278  static const ndarray<unsigned int, 6, 3> table = {{{{0, 2, 1}},
1279  {{0, 1, 2}},
1280  {{2, 1, 0}},
1281  {{1, 2, 0}},
1282  {{1, 0, 2}},
1283  {{2, 0, 1}}}};
1284 
1285  return table[face_orientation][vertex];
1286  }
1287  }
1288  else if (*this == ReferenceCells::Wedge)
1289  {
1290  if (face > 1) // One of the quadrilateral faces
1291  {
1293  vertex,
1294  Utilities::get_bit(face_orientation, 0),
1295  Utilities::get_bit(face_orientation, 2),
1296  Utilities::get_bit(face_orientation, 1));
1297  }
1298  else // One of the triangular faces
1299  {
1300  static const ndarray<unsigned int, 6, 3> table = {{{{0, 2, 1}},
1301  {{0, 1, 2}},
1302  {{2, 1, 0}},
1303  {{1, 2, 0}},
1304  {{1, 0, 2}},
1305  {{2, 0, 1}}}};
1306 
1307  return table[face_orientation][vertex];
1308  }
1309  }
1310  else if (*this == ReferenceCells::Hexahedron)
1311  {
1313  vertex,
1314  Utilities::get_bit(face_orientation, 0),
1315  Utilities::get_bit(face_orientation, 2),
1316  Utilities::get_bit(face_orientation, 1));
1317  }
1318 
1319  Assert(false, ExcNotImplemented());
1321 }
1322 
1323 
1324 
1325 inline unsigned int
1327  const unsigned int line,
1328  const unsigned int face,
1329  const unsigned char face_orientation) const
1330 {
1331  AssertIndexRange(face, n_faces());
1333 
1334  if (*this == ReferenceCells::Vertex)
1335  {
1336  Assert(false, ExcNotImplemented());
1337  }
1338  else if (*this == ReferenceCells::Line)
1339  {
1340  Assert(false, ExcNotImplemented());
1341  }
1342  else if (*this == ReferenceCells::Triangle)
1343  {
1344  Assert(false, ExcNotImplemented());
1345  }
1346  else if (*this == ReferenceCells::Quadrilateral)
1347  {
1348  Assert(false, ExcNotImplemented());
1349  }
1350  else if (*this == ReferenceCells::Tetrahedron)
1351  {
1352  static const ndarray<unsigned int, 6, 3> table = {{{{2, 1, 0}},
1353  {{0, 1, 2}},
1354  {{1, 0, 2}},
1355  {{1, 2, 0}},
1356  {{0, 2, 1}},
1357  {{2, 0, 1}}}};
1358 
1359  return table[face_orientation][line];
1360  }
1361  else if (*this == ReferenceCells::Pyramid)
1362  {
1363  if (face == 0) // The quadrilateral face
1364  {
1366  line,
1367  Utilities::get_bit(face_orientation, 0),
1368  Utilities::get_bit(face_orientation, 2),
1369  Utilities::get_bit(face_orientation, 1));
1370  }
1371  else // One of the triangular faces
1372  {
1373  static const ndarray<unsigned int, 6, 3> table = {{{{2, 1, 0}},
1374  {{0, 1, 2}},
1375  {{1, 0, 2}},
1376  {{1, 2, 0}},
1377  {{0, 2, 1}},
1378  {{2, 0, 1}}}};
1379 
1380  return table[face_orientation][line];
1381  }
1382  }
1383  else if (*this == ReferenceCells::Wedge)
1384  {
1385  if (face > 1) // One of the quadrilateral faces
1386  {
1388  line,
1389  Utilities::get_bit(face_orientation, 0),
1390  Utilities::get_bit(face_orientation, 2),
1391  Utilities::get_bit(face_orientation, 1));
1392  }
1393  else // One of the triangular faces
1394  {
1395  static const ndarray<unsigned int, 6, 3> table = {{{{2, 1, 0}},
1396  {{0, 1, 2}},
1397  {{1, 0, 2}},
1398  {{1, 2, 0}},
1399  {{0, 2, 1}},
1400  {{2, 0, 1}}}};
1401 
1402  return table[face_orientation][line];
1403  }
1404  }
1405  else if (*this == ReferenceCells::Hexahedron)
1406  {
1408  line,
1409  Utilities::get_bit(face_orientation, 0),
1410  Utilities::get_bit(face_orientation, 2),
1411  Utilities::get_bit(face_orientation, 1));
1412  }
1413 
1414  Assert(false, ExcNotImplemented());
1416 }
1417 
1418 
1419 
1420 namespace ReferenceCells
1421 {
1422  template <int dim>
1423  inline constexpr const ReferenceCell &
1425  {
1426  switch (dim)
1427  {
1428  case 0:
1429  return ReferenceCells::Vertex;
1430  case 1:
1431  return ReferenceCells::Line;
1432  case 2:
1433  return ReferenceCells::Triangle;
1434  case 3:
1436  default:
1437  Assert(false, ExcNotImplemented());
1438  return ReferenceCells::Invalid;
1439  }
1440  }
1441 
1442 
1443 
1444  template <int dim>
1445  inline constexpr const ReferenceCell &
1447  {
1448  switch (dim)
1449  {
1450  case 0:
1451  return ReferenceCells::Vertex;
1452  case 1:
1453  return ReferenceCells::Line;
1454  case 2:
1456  case 3:
1458  default:
1459  Assert(false, ExcNotImplemented());
1460  return ReferenceCells::Invalid;
1461  }
1462  }
1463 } // namespace ReferenceCells
1464 
1465 
1466 inline ReferenceCell
1467 ReferenceCell::n_vertices_to_type(const int dim, const unsigned int n_vertices)
1468 {
1469  AssertIndexRange(dim, 4);
1470  AssertIndexRange(n_vertices, 9);
1471 
1472  const auto X = ReferenceCells::Invalid;
1473  static const ndarray<ReferenceCell, 4, 9> table = {
1474  {// dim 0
1475  {{X, ReferenceCells::Vertex, X, X, X, X, X, X, X}},
1476  // dim 1
1477  {{X, X, ReferenceCells::Line, X, X, X, X, X, X}},
1478  // dim 2
1479  {{X,
1480  X,
1481  X,
1484  X,
1485  X,
1486  X,
1487  X}},
1488  // dim 3
1489  {{X,
1490  X,
1491  X,
1492  X,
1496  X,
1498  Assert(table[dim][n_vertices] != ReferenceCells::Invalid,
1499  ExcMessage("The combination of dim = " + std::to_string(dim) +
1500  " and n_vertices = " + std::to_string(n_vertices) +
1501  " does not correspond to a known reference cell type."));
1502  return table[dim][n_vertices];
1503 }
1504 
1505 
1506 
1507 template <int dim>
1508 inline double
1510  const unsigned int i) const
1511 {
1513  if (*this == ReferenceCells::get_hypercube<dim>())
1515 
1516  if (*this ==
1517  ReferenceCells::Triangle) // see also
1518  // BarycentricPolynomials<2>::compute_value
1519  {
1520  switch (i)
1521  {
1522  case 0:
1523  return 1.0 - xi[std::min(0, dim - 1)] - xi[std::min(1, dim - 1)];
1524  case 1:
1525  return xi[std::min(0, dim - 1)];
1526  case 2:
1527  return xi[std::min(1, dim - 1)];
1528  }
1529  }
1530 
1531  if (*this ==
1532  ReferenceCells::Tetrahedron) // see also
1533  // BarycentricPolynomials<3>::compute_value
1534  {
1535  switch (i)
1536  {
1537  case 0:
1538  return 1.0 - xi[std::min(0, dim - 1)] - xi[std::min(1, dim - 1)] -
1539  xi[std::min(2, dim - 1)];
1540  case 1:
1541  return xi[std::min(0, dim - 1)];
1542  case 2:
1543  return xi[std::min(1, dim - 1)];
1544  case 3:
1545  return xi[std::min(2, dim - 1)];
1546  }
1547  }
1548 
1549  if (*this ==
1550  ReferenceCells::Wedge) // see also
1551  // ScalarLagrangePolynomialWedge::compute_value
1552  {
1554  .d_linear_shape_function<2>(Point<2>(xi[std::min(0, dim - 1)],
1555  xi[std::min(1, dim - 1)]),
1556  i % 3) *
1558  .d_linear_shape_function<1>(Point<1>(xi[std::min(2, dim - 1)]),
1559  i / 3);
1560  }
1561 
1562  if (*this ==
1563  ReferenceCells::Pyramid) // see also
1564  // ScalarLagrangePolynomialPyramid::compute_value
1565  {
1566  const double Q14 = 0.25;
1567  double ration;
1568 
1569  const double r = xi[std::min(0, dim - 1)];
1570  const double s = xi[std::min(1, dim - 1)];
1571  const double t = xi[std::min(2, dim - 1)];
1572 
1573  if (fabs(t - 1.0) > 1.0e-14)
1574  {
1575  ration = (r * s * t) / (1.0 - t);
1576  }
1577  else
1578  {
1579  ration = 0.0;
1580  }
1581 
1582  if (i == 0)
1583  return Q14 * ((1.0 - r) * (1.0 - s) - t + ration);
1584  if (i == 1)
1585  return Q14 * ((1.0 + r) * (1.0 - s) - t - ration);
1586  if (i == 2)
1587  return Q14 * ((1.0 - r) * (1.0 + s) - t - ration);
1588  if (i == 3)
1589  return Q14 * ((1.0 + r) * (1.0 + s) - t + ration);
1590  else
1591  return t;
1592  }
1593 
1594  Assert(false, ExcNotImplemented());
1595 
1596  return 0.0;
1597 }
1598 
1599 
1600 
1601 template <int dim>
1602 inline Tensor<1, dim>
1604  const unsigned int i) const
1605 {
1607  if (*this == ReferenceCells::get_hypercube<dim>())
1609 
1610  if (*this ==
1611  ReferenceCells::Triangle) // see also
1612  // BarycentricPolynomials<2>::compute_grad
1613  {
1614  switch (i)
1615  {
1616  case 0:
1617  return Point<dim>(-1.0, -1.0);
1618  case 1:
1619  return Point<dim>(+1.0, +0.0);
1620  case 2:
1621  return Point<dim>(+0.0, +1.0);
1622  }
1623  }
1624 
1625  Assert(false, ExcNotImplemented());
1626 
1627  return Point<dim>(+0.0, +0.0, +0.0);
1628 }
1629 
1630 
1631 template <int dim>
1632 inline Tensor<1, dim>
1633 ReferenceCell::unit_tangential_vectors(const unsigned int face_no,
1634  const unsigned int i) const
1635 {
1637  AssertIndexRange(i, dim - 1);
1638 
1639  if (*this == ReferenceCells::get_hypercube<dim>())
1640  {
1642  return GeometryInfo<dim>::unit_tangential_vectors[face_no][i];
1643  }
1644  else if (*this == ReferenceCells::Triangle)
1645  {
1646  AssertIndexRange(face_no, 3);
1647  static const std::array<Tensor<1, dim>, 3> table = {
1648  {Point<dim>(1, 0),
1649  Point<dim>(-std::sqrt(0.5), +std::sqrt(0.5)),
1650  Point<dim>(0, -1)}};
1651 
1652  return table[face_no];
1653  }
1654  else if (*this == ReferenceCells::Tetrahedron)
1655  {
1656  AssertIndexRange(face_no, 4);
1657  static const ndarray<Tensor<1, dim>, 4, 2> table = {
1658  {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
1659  {{Point<dim>(1, 0, 0), Point<dim>(0, 0, 1)}},
1660  {{Point<dim>(0, 0, 1), Point<dim>(0, 1, 0)}},
1661  {{Point<dim>(-std::pow(1.0 / 3.0, 1.0 / 4.0),
1662  +std::pow(1.0 / 3.0, 1.0 / 4.0),
1663  0),
1664  Point<dim>(-std::pow(1.0 / 3.0, 1.0 / 4.0),
1665  0,
1666  +std::pow(1.0 / 3.0, 1.0 / 4.0))}}}};
1667 
1668  return table[face_no][i];
1669  }
1670  else if (*this == ReferenceCells::Wedge)
1671  {
1672  AssertIndexRange(face_no, 5);
1673  static const ndarray<Tensor<1, dim>, 5, 2> table = {
1674  {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
1675  {{Point<dim>(1, 0, 0), Point<dim>(0, 1, 0)}},
1676  {{Point<dim>(1, 0, 0), Point<dim>(0, 0, 1)}},
1677  {{Point<dim>(-1 / std::sqrt(2.0), +1 / std::sqrt(2.0), 0),
1678  Point<dim>(0, 0, 1)}},
1679  {{Point<dim>(0, 0, 1), Point<dim>(0, 1, 0)}}}};
1680 
1681  return table[face_no][i];
1682  }
1683  else if (*this == ReferenceCells::Pyramid)
1684  {
1685  AssertIndexRange(face_no, 5);
1686  static const ndarray<Tensor<1, dim>, 5, 2> table = {
1687  {{{Point<dim>(0, 1, 0), Point<dim>(1, 0, 0)}},
1688  {{Point<dim>(+1.0 / sqrt(2.0), 0, +1.0 / sqrt(2.0)),
1689  Point<dim>(0, 1, 0)}},
1690  {{Point<dim>(+1.0 / sqrt(2.0), 0, -1.0 / sqrt(2.0)),
1691  Point<dim>(0, 1, 0)}},
1692  {{Point<dim>(1, 0, 0),
1693  Point<dim>(0, +1.0 / sqrt(2.0), +1.0 / sqrt(2.0))}},
1694  {{Point<dim>(1, 0, 0),
1695  Point<dim>(0, +1.0 / sqrt(2.0), -1.0 / sqrt(2.0))}}}};
1696 
1697  return table[face_no][i];
1698  }
1699 
1700  Assert(false, ExcNotImplemented());
1701 
1702  return {};
1703 }
1704 
1705 
1706 
1707 template <int dim>
1708 inline Tensor<1, dim>
1709 ReferenceCell::unit_normal_vectors(const unsigned int face_no) const
1710 {
1711  AssertDimension(dim, this->get_dimension());
1712 
1713  if (is_hyper_cube())
1714  {
1716  return GeometryInfo<dim>::unit_normal_vector[face_no];
1717  }
1718  else if (dim == 2)
1719  {
1721 
1722  // Return the rotated vector
1723  return cross_product_2d(unit_tangential_vectors<dim>(face_no, 0));
1724  }
1725  else if (dim == 3)
1726  {
1727  return cross_product_3d(unit_tangential_vectors<dim>(face_no, 0),
1728  unit_tangential_vectors<dim>(face_no, 1));
1729  }
1730 
1731  Assert(false, ExcNotImplemented());
1732 
1733  return {};
1734 }
1735 
1736 
1737 
1738 inline bool
1740  const unsigned int line,
1741  const unsigned char face_orientation_raw,
1742  const unsigned char line_orientation) const
1743 {
1744  if (*this == ReferenceCells::Hexahedron)
1745  {
1746  static const bool bool_table[2][2][2][2] = {
1747  {{{true, false}, // lines 0/1, face_orientation=false,
1748  // face_flip=false, face_rotation=false and true
1749  {false, true}}, // lines 0/1, face_orientation=false,
1750  // face_flip=true, face_rotation=false and true
1751  {{true, true}, // lines 0/1, face_orientation=true,
1752  // face_flip=false, face_rotation=false and true
1753  {false, false}}}, // lines 0/1, face_orientation=true,
1754  // face_flip=true, face_rotation=false and true
1755 
1756  {{{true, true}, // lines 2/3 ...
1757  {false, false}},
1758  {{true, false}, {false, true}}}};
1759 
1760  const bool face_orientation = Utilities::get_bit(face_orientation_raw, 0);
1761  const bool face_flip = Utilities::get_bit(face_orientation_raw, 2);
1762  const bool face_rotation = Utilities::get_bit(face_orientation_raw, 1);
1763 
1764  return (static_cast<bool>(line_orientation) ==
1765  bool_table[line / 2][face_orientation][face_flip][face_rotation]);
1766  }
1767  else
1768  // TODO: This might actually be wrong for some of the other
1769  // kinds of objects. We should check this
1770  return true;
1771 }
1772 
1773 
1774 
1775 namespace internal
1776 {
1777  template <typename T, std::size_t N>
1779  {
1780  public:
1784  NoPermutation(const ::ReferenceCell &entity_type,
1785  const std::array<T, N> & vertices_0,
1786  const std::array<T, N> & vertices_1)
1787  : entity_type(entity_type)
1788  , vertices_0(vertices_0)
1789  , vertices_1(vertices_1)
1790  {}
1791 
1795  virtual ~NoPermutation() noexcept override = default;
1796 
1800  virtual void
1801  print_info(std::ostream &out) const override
1802  {
1803  out << "[";
1804 
1805  const unsigned int n_vertices = entity_type.n_vertices();
1806 
1807  for (unsigned int i = 0; i < n_vertices; ++i)
1808  {
1809  out << vertices_0[i];
1810  if (i + 1 != n_vertices)
1811  out << ",";
1812  }
1813 
1814  out << "] is not a permutation of [";
1815 
1816  for (unsigned int i = 0; i < n_vertices; ++i)
1817  {
1818  out << vertices_1[i];
1819  if (i + 1 != n_vertices)
1820  out << ",";
1821  }
1822 
1823  out << "]." << std::endl;
1824  }
1825 
1829  const ::ReferenceCell entity_type;
1830 
1834  const std::array<T, N> vertices_0;
1835 
1839  const std::array<T, N> vertices_1;
1840  };
1841 } // namespace internal
1842 
1843 
1844 
1845 template <typename T, std::size_t N>
1846 inline unsigned char
1847 ReferenceCell::compute_orientation(const std::array<T, N> &vertices_0,
1848  const std::array<T, N> &vertices_1) const
1849 {
1850  AssertIndexRange(n_vertices(), N + 1);
1851  if (*this == ReferenceCells::Line)
1852  {
1853  const std::array<T, 2> i{{vertices_0[0], vertices_0[1]}};
1854  const std::array<T, 2> j{{vertices_1[0], vertices_1[1]}};
1855 
1856  // line_orientation=true
1857  if (i == std::array<T, 2>{{j[0], j[1]}})
1858  return 1;
1859 
1860  // line_orientation=false
1861  if (i == std::array<T, 2>{{j[1], j[0]}})
1862  return 0;
1863  }
1864  else if (*this == ReferenceCells::Triangle)
1865  {
1866  const std::array<T, 3> i{{vertices_0[0], vertices_0[1], vertices_0[2]}};
1867  const std::array<T, 3> j{{vertices_1[0], vertices_1[1], vertices_1[2]}};
1868 
1869  // face_orientation=true, face_rotation=false, face_flip=false
1870  if (i == std::array<T, 3>{{j[0], j[1], j[2]}})
1871  return 1;
1872 
1873  // face_orientation=true, face_rotation=true, face_flip=false
1874  if (i == std::array<T, 3>{{j[1], j[2], j[0]}})
1875  return 3;
1876 
1877  // face_orientation=true, face_rotation=false, face_flip=true
1878  if (i == std::array<T, 3>{{j[2], j[0], j[1]}})
1879  return 5;
1880 
1881  // face_orientation=false, face_rotation=false, face_flip=false
1882  if (i == std::array<T, 3>{{j[0], j[2], j[1]}})
1883  return 0;
1884 
1885  // face_orientation=false, face_rotation=true, face_flip=false
1886  if (i == std::array<T, 3>{{j[2], j[1], j[0]}})
1887  return 2;
1888 
1889  // face_orientation=false, face_rotation=false, face_flip=true
1890  if (i == std::array<T, 3>{{j[1], j[0], j[2]}})
1891  return 4;
1892  }
1893  else if (*this == ReferenceCells::Quadrilateral)
1894  {
1895  const std::array<T, 4> i{
1896  {vertices_0[0], vertices_0[1], vertices_0[2], vertices_0[3]}};
1897  const std::array<T, 4> j{
1898  {vertices_1[0], vertices_1[1], vertices_1[2], vertices_1[3]}};
1899 
1900  // face_orientation=true, face_rotation=false, face_flip=false
1901  if (i == std::array<T, 4>{{j[0], j[1], j[2], j[3]}})
1902  return 1;
1903 
1904  // face_orientation=true, face_rotation=true, face_flip=false
1905  if (i == std::array<T, 4>{{j[2], j[0], j[3], j[1]}})
1906  return 3;
1907 
1908  // face_orientation=true, face_rotation=false, face_flip=true
1909  if (i == std::array<T, 4>{{j[3], j[2], j[1], j[0]}})
1910  return 5;
1911 
1912  // face_orientation=true, face_rotation=true, face_flip=true
1913  if (i == std::array<T, 4>{{j[1], j[3], j[0], j[2]}})
1914  return 7;
1915 
1916  // face_orientation=false, face_rotation=false, face_flip=false
1917  if (i == std::array<T, 4>{{j[0], j[2], j[1], j[3]}})
1918  return 0;
1919 
1920  // face_orientation=false, face_rotation=true, face_flip=false
1921  if (i == std::array<T, 4>{{j[2], j[3], j[0], j[1]}})
1922  return 2;
1923 
1924  // face_orientation=false, face_rotation=false, face_flip=true
1925  if (i == std::array<T, 4>{{j[3], j[1], j[2], j[0]}})
1926  return 4;
1927 
1928  // face_orientation=false, face_rotation=true, face_flip=true
1929  if (i == std::array<T, 4>{{j[1], j[0], j[3], j[2]}})
1930  return 6;
1931  }
1932 
1933  Assert(false, (internal::NoPermutation<T, N>(*this, vertices_0, vertices_1)));
1934 
1935  return -1;
1936 }
1937 
1938 
1939 
1940 template <typename T, std::size_t N>
1941 inline std::array<T, N>
1943  const std::array<T, N> &vertices,
1944  const unsigned int orientation) const
1945 {
1946  std::array<T, 4> temp;
1947 
1948  if (*this == ReferenceCells::Line)
1949  {
1950  switch (orientation)
1951  {
1952  case 1:
1953  temp = {{vertices[0], vertices[1]}};
1954  break;
1955  case 0:
1956  temp = {{vertices[1], vertices[0]}};
1957  break;
1958  default:
1959  Assert(false, ExcNotImplemented());
1960  }
1961  }
1962  else if (*this == ReferenceCells::Triangle)
1963  {
1964  switch (orientation)
1965  {
1966  case 1:
1967  temp = {{vertices[0], vertices[1], vertices[2]}};
1968  break;
1969  case 3:
1970  temp = {{vertices[1], vertices[2], vertices[0]}};
1971  break;
1972  case 5:
1973  temp = {{vertices[2], vertices[0], vertices[1]}};
1974  break;
1975  case 0:
1976  temp = {{vertices[0], vertices[2], vertices[1]}};
1977  break;
1978  case 2:
1979  temp = {{vertices[2], vertices[1], vertices[0]}};
1980  break;
1981  case 4:
1982  temp = {{vertices[1], vertices[0], vertices[2]}};
1983  break;
1984  default:
1985  Assert(false, ExcNotImplemented());
1986  }
1987  }
1988  else if (*this == ReferenceCells::Quadrilateral)
1989  {
1990  switch (orientation)
1991  {
1992  case 1:
1993  temp = {{vertices[0], vertices[1], vertices[2], vertices[3]}};
1994  break;
1995  case 3:
1996  temp = {{vertices[2], vertices[0], vertices[3], vertices[1]}};
1997  break;
1998  case 5:
1999  temp = {{vertices[3], vertices[2], vertices[1], vertices[0]}};
2000  break;
2001  case 7:
2002  temp = {{vertices[1], vertices[3], vertices[0], vertices[2]}};
2003  break;
2004  case 0:
2005  temp = {{vertices[0], vertices[2], vertices[1], vertices[3]}};
2006  break;
2007  case 2:
2008  temp = {{vertices[2], vertices[3], vertices[0], vertices[1]}};
2009  break;
2010  case 4:
2011  temp = {{vertices[3], vertices[1], vertices[2], vertices[0]}};
2012  break;
2013  case 6:
2014  temp = {{vertices[1], vertices[0], vertices[3], vertices[2]}};
2015  break;
2016  default:
2017  Assert(false, ExcNotImplemented());
2018  }
2019  }
2020  else
2021  {
2022  AssertThrow(false, ExcNotImplemented());
2023  }
2024 
2025  std::array<T, N> temp_;
2026  std::copy_n(temp.begin(), N, temp_.begin());
2027 
2028  return temp_;
2029 }
2030 
2031 
2033 
2034 #endif
unsigned int standard_to_real_face_line(const unsigned int line, const unsigned int face, const unsigned char face_orientation) const
constexpr const ReferenceCell Pyramid
static const unsigned int invalid_unsigned_int
Definition: types.h:196
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)
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1623
bool is_hyper_cube() const
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)
ArrayView< const unsigned int > faces_for_given_vertex(const unsigned int vertex_index) const
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)
constexpr const ReferenceCell & get_hypercube()
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char face_orientation) const
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
const Mapping< dim, spacedim > & get_default_linear_mapping(const Triangulation< dim, spacedim > &triangulation)
Definition: mapping.cc:240
bool is_simplex() const
#define AssertIndexRange(index, range)
Definition: exceptions.h:1691
constexpr const ReferenceCell Triangle
static ReferenceCell n_vertices_to_type(const int dim, const unsigned int n_vertices)
constexpr ::ReferenceCell make_reference_cell_from_int(const std::uint8_t kind)
static const char U
Tensor< 1, dim > unit_normal_vectors(const unsigned int face_no) const
static Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i)
double d_linear_shape_function(const Point< dim > &xi, const unsigned int i) const
#define AssertThrow(cond, exc)
Definition: exceptions.h:1576
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
constexpr Tensor< 1, dim, typename ProductType< Number1, Number2 >::type > cross_product_3d(const Tensor< 1, dim, Number1 > &src1, const Tensor< 1, dim, Number2 > &src2)
Definition: tensor.h:2443
void serialize(Archive &archive, const unsigned int)
constexpr bool operator==(const ReferenceCell &type) const
constexpr const ReferenceCell Line
static std::array< unsigned int, 2 > standard_quad_vertex_to_line_vertex_index(const unsigned int vertex)
constexpr const ReferenceCell Wedge
Tensor< 1, dim > unit_tangential_vectors(const unsigned int face_no, const unsigned int i) const
static ::ExceptionBase & ExcMessage(std::string arg1)
NoPermutation(const ::ReferenceCell &entity_type, const std::array< T, N > &vertices_0, const std::array< T, N > &vertices_1)
bool standard_vs_true_line_orientation(const unsigned int line, const unsigned char face_orientation, const unsigned char line_orientation) const
constexpr const ReferenceCell Tetrahedron
#define Assert(cond, exc)
Definition: exceptions.h:1466
std::uint8_t kind
constexpr bool operator!=(const ReferenceCell &type) const
Abstract base class for mapping classes.
Definition: mapping.h:303
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
std::array< unsigned int, 2 > standard_line_to_face_and_line_index(const unsigned int line) const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:396
std::array< T, N > permute_according_orientation(const std::array< T, N > &vertices, const unsigned int orientation) const
static unsigned int standard_to_real_line_vertex(const unsigned int vertex, const bool line_orientation=true)
std::string to_string(const T &t)
Definition: patterns.h:2329
constexpr const ReferenceCell Hexahedron
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
Point< 3 > vertices[4]
unsigned int n_vertices() const
unsigned int child_cell_on_face(const unsigned int face_n, const unsigned int subface_n) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
static unsigned int standard_to_real_face_vertex(const unsigned int vertex, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
Expression fabs(const Expression &x)
unsigned char compute_orientation(const std::array< T, N > &vertices_0, const std::array< T, N > &vertices_1) const
virtual void print_info(std::ostream &out) const override
unsigned int get_dimension() const
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const
unsigned int n_lines() const
unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const unsigned char face_orientation) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
static double d_linear_shape_function(const Point< dim > &xi, const unsigned int i)
constexpr const ReferenceCell Invalid
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:395
T min(const T &t, const MPI_Comm &mpi_communicator)
constexpr const ReferenceCell Quadrilateral
static const char N
constexpr const ReferenceCell Vertex
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const
static std::array< unsigned int, 2 > standard_hex_vertex_to_quad_vertex_index(const unsigned int vertex)
static ::ExceptionBase & ExcNotImplemented()
Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i) const
unsigned int n_faces() const
constexpr ReferenceCell()
const std::array< T, N > vertices_1
const ::ReferenceCell entity_type
constexpr Tensor< 1, dim, Number > cross_product_2d(const Tensor< 1, dim, Number > &src)
Definition: tensor.h:2418
constexpr const ReferenceCell & get_simplex()
ReferenceCell face_reference_cell(const unsigned int face_no) const
static std::array< unsigned int, 2 > standard_hex_line_to_quad_line_index(const unsigned int line)
const std::array< T, N > vertices_0
typename internal::ndarray::HelperArray< T, Ns... >::type ndarray
Definition: ndarray.h:108
#define DEAL_II_CONSTEXPR
Definition: config.h:168
bool get_bit(const unsigned char number, const unsigned int n)
Definition: utilities.h:1376
static ::ExceptionBase & ExcInternalError()