Reference documentation for deal.II version Git e12e976e47 2020-11-24 09:52:11 +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 
22 
23 
25 
26 
30 namespace ReferenceCell
31 {
35  enum class Type : std::uint8_t
36  {
37  Vertex = 0,
38  Line = 1,
39  Tri = 2,
40  Quad = 3,
41  Tet = 4,
42  Pyramid = 5,
43  Wedge = 6,
44  Hex = 7,
45  Invalid = static_cast<std::uint8_t>(-1)
46  };
47 
52  inline Type
53  get_simplex(const unsigned int dim)
54  {
55  switch (dim)
56  {
57  case 0:
58  return Type::Vertex;
59  case 1:
60  return Type::Line;
61  case 2:
62  return Type::Tri;
63  case 3:
64  return Type::Tet;
65  default:
66  Assert(false, ExcNotImplemented());
67  return Type::Invalid;
68  }
69  }
70 
75  inline Type
76  get_hypercube(const unsigned int dim)
77  {
78  switch (dim)
79  {
80  case 0:
81  return Type::Vertex;
82  case 1:
83  return Type::Line;
84  case 2:
85  return Type::Quad;
86  case 3:
87  return Type::Hex;
88  default:
89  Assert(false, ExcNotImplemented());
90  return Type::Invalid;
91  }
92  }
93 
98  inline Type
99  n_vertices_to_type(const int dim, const unsigned int n_vertices)
100  {
101  AssertIndexRange(dim, 4);
102  AssertIndexRange(n_vertices, 9);
103  const auto X = Type::Invalid;
104 
105  static constexpr std::array<std::array<ReferenceCell::Type, 9>, 4> table = {
106  {// dim 0
107  {{X, Type::Vertex, X, X, X, X, X, X, X}},
108  // dim 1
109  {{X, X, Type::Line, X, X, X, X, X, X}},
110  // dim 2
111  {{X, X, X, Type::Tri, Type::Quad, X, X, X, X}},
112  // dim 3
113  {{X, X, X, X, Type::Tet, Type::Pyramid, Type::Wedge, X, Type::Hex}}}};
114  Assert(table[dim][n_vertices] != Type::Invalid,
115  ExcMessage("The combination of dim = " + std::to_string(dim) +
116  " and n_vertices = " + std::to_string(n_vertices) +
117  " does not correspond to a known reference cell type."));
118  return table[dim][n_vertices];
119  }
120 
121  namespace internal
122  {
126  inline static bool
127  get_bit(const unsigned char number, const unsigned int n)
128  {
129  AssertIndexRange(n, 8);
130 
131  // source:
132  // https://stackoverflow.com/questions/47981/how-do-you-set-clear-and-toggle-a-single-bit
133  // "Checking a bit"
134  return (number >> n) & 1U;
135  }
136 
137 
138 
142  inline static void
143  set_bit(unsigned char &number, const unsigned int n, const bool x)
144  {
145  AssertIndexRange(n, 8);
146 
147  // source:
148  // https://stackoverflow.com/questions/47981/how-do-you-set-clear-and-toggle-a-single-bit
149  // "Changing the nth bit to x"
150  number ^= (-static_cast<unsigned char>(x) ^ number) & (1U << n);
151  }
152 
156  namespace Info
157  {
163  struct Base
164  {
168  virtual ~Base() = default;
169 
173  virtual unsigned int
174  n_vertices() const
175  {
176  Assert(false, ExcNotImplemented());
177  return 0;
178  }
179 
183  virtual unsigned int
184  n_lines() const
185  {
186  Assert(false, ExcNotImplemented());
187  return 0;
188  }
189 
190 
194  virtual unsigned int
195  n_faces() const
196  {
197  Assert(false, ExcNotImplemented());
198  return 0;
199  }
200 
207  {
208  return {0U, n_vertices()};
209  }
210 
216  line_indices() const
217  {
218  return {0U, n_lines()};
219  }
220 
226  face_indices() const
227  {
228  return {0U, n_faces()};
229  }
230 
235  virtual std::array<unsigned int, 2>
237  const unsigned int vertex) const
238  {
239  Assert(false, ExcNotImplemented());
240 
241  (void)vertex;
242 
243  return {{0u, 0u}};
244  }
245 
249  virtual std::array<unsigned int, 2>
250  standard_line_to_face_and_line_index(const unsigned int line) const
251  {
252  Assert(false, ExcNotImplemented());
253 
254  (void)line;
255 
256  return {{0, 0}};
257  }
258 
262  virtual unsigned int
263  standard_to_real_face_vertex(const unsigned int vertex,
264  const unsigned int face,
265  const unsigned char face_orientation) const
266  {
267  Assert(false, ExcNotImplemented());
268 
269  (void)vertex;
270  (void)face;
271  (void)face_orientation;
272 
273  return 0;
274  }
275 
279  virtual unsigned int
280  standard_to_real_face_line(const unsigned int line,
281  const unsigned int face,
282  const unsigned char face_orientation) const
283  {
284  Assert(false, ExcNotImplemented());
285 
286  (void)line;
287  (void)face;
288  (void)face_orientation;
289 
290  return 0;
291  }
292 
296  virtual bool
298  const unsigned int line,
299  const unsigned char face_orientation,
300  const unsigned char line_orientation) const
301  {
302  Assert(false, ExcNotImplemented());
303 
304  (void)line;
305  (void)face_orientation;
306  (void)line_orientation;
307 
308  return true;
309  }
310 
314  virtual ReferenceCell::Type
315  face_reference_cell_type(const unsigned int face_no) const
316  {
317  Assert(false, ExcNotImplemented());
318  (void)face_no;
319 
321  }
322 
326  virtual unsigned int
327  face_to_cell_lines(const unsigned int face,
328  const unsigned int line,
329  const unsigned char face_orientation) const
330  {
331  Assert(false, ExcNotImplemented());
332  (void)face;
333  (void)line;
334  (void)face_orientation;
335 
336  return 0;
337  }
338 
342  virtual unsigned int
343  face_to_cell_vertices(const unsigned int face,
344  const unsigned int vertex,
345  const unsigned char face_orientation) const
346  {
347  Assert(false, ExcNotImplemented());
348  (void)face;
349  (void)vertex;
350  (void)face_orientation;
351 
352  return 0;
353  }
354 
358  virtual unsigned int
359  exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const
360  {
361  Assert(false, ExcNotImplemented());
362  (void)vertex_n;
363 
364  return 0;
365  }
366 
370  virtual unsigned int
371  exodusii_face_to_deal_face(const unsigned int face_n) const
372  {
373  Assert(false, ExcNotImplemented());
374  (void)face_n;
375 
376  return 0;
377  }
378  };
379 
380 
384  template <int dim>
386  {
387  unsigned int
388  n_vertices() const override
389  {
391  }
392 
393  unsigned int
394  n_lines() const override
395  {
397  }
398 
399  unsigned int
400  n_faces() const override
401  {
403  }
404 
405  unsigned int
406  face_to_cell_lines(const unsigned int face,
407  const unsigned int line,
408  const unsigned char face_orientation) const override
409  {
411  face,
412  line,
413  get_bit(face_orientation, 0),
414  get_bit(face_orientation, 2),
415  get_bit(face_orientation, 1));
416  }
417 
418  unsigned int
420  const unsigned int face,
421  const unsigned int vertex,
422  const unsigned char face_orientation) const override
423  {
425  face,
426  vertex,
427  get_bit(face_orientation, 0),
428  get_bit(face_orientation, 2),
429  get_bit(face_orientation, 1));
430  }
431  };
432 
433 
434 
435  /*
436  * Vertex.
437  */
438  struct Vertex : public TensorProductBase<0>
439  {
441  face_reference_cell_type(const unsigned int face_no) const override
442  {
443  (void)face_no;
445  }
446 
447  virtual unsigned int
448  exodusii_face_to_deal_face(const unsigned int face_n) const override
449  {
450  (void)face_n;
451  AssertIndexRange(face_n, n_faces());
452 
453  return 0;
454  }
455  };
456 
457 
458 
459  /*
460  * Line.
461  */
462  struct Line : public TensorProductBase<1>
463  {
465  face_reference_cell_type(const unsigned int face_no) const override
466  {
467  (void)face_no;
469  }
470 
471  virtual unsigned int
473  const unsigned int vertex_n) const override
474  {
475  AssertIndexRange(vertex_n, n_vertices());
476  return vertex_n;
477  }
478 
479  virtual unsigned int
480  exodusii_face_to_deal_face(const unsigned int face_n) const override
481  {
482  AssertIndexRange(face_n, n_faces());
483  return face_n;
484  }
485  };
486 
487 
488 
492  struct Tri : public Base
493  {
494  unsigned int
495  n_vertices() const override
496  {
497  return 3;
498  }
499 
500  unsigned int
501  n_lines() const override
502  {
503  return 3;
504  }
505 
506  unsigned int
507  n_faces() const override
508  {
509  return this->n_lines();
510  }
511 
512  std::array<unsigned int, 2>
514  const unsigned int vertex) const override
515  {
516  AssertIndexRange(vertex, 3);
517 
518  static const std::array<std::array<unsigned int, 2>, 3> table = {
519  {{{0, 0}}, {{0, 1}}, {{1, 1}}}};
520 
521  return table[vertex];
522  }
523 
524  unsigned int
526  const unsigned int vertex,
527  const unsigned int face,
528  const unsigned char line_orientation) const override
529  {
530  (void)face;
531 
532  static const std::array<std::array<unsigned int, 2>, 2> table = {
533  {{{1, 0}}, {{0, 1}}}};
534 
535  return table[line_orientation][vertex];
536  }
537 
539  face_reference_cell_type(const unsigned int face_no) const override
540  {
541  (void)face_no;
542 
543  AssertIndexRange(face_no, n_faces());
544 
546  }
547 
548  unsigned int
549  face_to_cell_lines(const unsigned int face,
550  const unsigned int line,
551  const unsigned char face_orientation) const override
552  {
553  AssertIndexRange(face, n_faces());
554  AssertDimension(line, 0);
555 
556  (void)line;
557  (void)face_orientation;
558 
559  return face;
560  }
561 
562  unsigned int
564  const unsigned int face,
565  const unsigned int vertex,
566  const unsigned char face_orientation) const override
567  {
568  static const std::array<std::array<unsigned int, 2>, 3> table = {
569  {{{0, 1}}, {{1, 2}}, {{2, 0}}}};
570 
571  return table[face][face_orientation ? vertex : (1 - vertex)];
572  }
573 
574  virtual unsigned int
576  const unsigned int vertex_n) const override
577  {
578  AssertIndexRange(vertex_n, n_vertices());
579  return vertex_n;
580  }
581 
582  virtual unsigned int
583  exodusii_face_to_deal_face(const unsigned int face_n) const override
584  {
585  AssertIndexRange(face_n, n_faces());
586  return face_n;
587  }
588  };
589 
590 
591 
595  struct Quad : public TensorProductBase<2>
596  {
597  std::array<unsigned int, 2>
599  const unsigned int vertex) const override
600  {
602  vertex);
603  }
604 
605  unsigned int
607  const unsigned int vertex,
608  const unsigned int face,
609  const unsigned char line_orientation) const override
610  {
611  (void)face;
612 
614  vertex, line_orientation);
615  }
616 
618  face_reference_cell_type(const unsigned int face_no) const override
619  {
620  (void)face_no;
622  }
623 
624  virtual unsigned int
626  const unsigned int vertex_n) const override
627  {
628  AssertIndexRange(vertex_n, n_vertices());
629  constexpr std::array<unsigned int, 4> exodus_to_deal{{0, 1, 3, 2}};
630  return exodus_to_deal[vertex_n];
631  }
632 
633  virtual unsigned int
634  exodusii_face_to_deal_face(const unsigned int face_n) const override
635  {
636  AssertIndexRange(face_n, n_faces());
637  constexpr std::array<unsigned int, 4> exodus_to_deal{{2, 1, 3, 0}};
638  return exodus_to_deal[face_n];
639  }
640  };
641 
642 
643 
647  struct Tet : public TensorProductBase<3>
648  {
649  unsigned int
650  n_vertices() const override
651  {
652  return 4;
653  }
654 
655  unsigned int
656  n_lines() const override
657  {
658  return 6;
659  }
660 
661  unsigned int
662  n_faces() const override
663  {
664  return 4;
665  }
666 
667  std::array<unsigned int, 2>
669  const unsigned int line) const override
670  {
671  static const std::array<unsigned int, 2> table[6] = {
672  {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}};
673 
674  return table[line];
675  }
676 
677  unsigned int
679  const unsigned int line,
680  const unsigned int face,
681  const unsigned char face_orientation) const override
682  {
683  (void)face;
684 
685  static const std::array<std::array<unsigned int, 3>, 6> table = {
686  {{{2, 1, 0}},
687  {{0, 1, 2}},
688  {{1, 2, 0}},
689  {{0, 2, 1}},
690  {{1, 0, 2}},
691  {{2, 0, 1}}}};
692 
693  return table[face_orientation][line];
694  }
695 
696  bool
698  const unsigned int line,
699  const unsigned char face_orientation_raw,
700  const unsigned char line_orientation) const override
701  {
702  (void)line;
703  (void)face_orientation_raw;
704 
705  return line_orientation;
706  }
707 
708  std::array<unsigned int, 2>
710  const unsigned int vertex) const override
711  {
712  AssertIndexRange(vertex, 4);
713 
714  static const std::array<unsigned int, 2> table[4] = {{{0, 0}},
715  {{0, 1}},
716  {{0, 2}},
717  {{1, 2}}};
718 
719  return table[vertex];
720  }
721 
722  unsigned int
724  const unsigned int vertex,
725  const unsigned int face,
726  const unsigned char face_orientation) const override
727  {
728  AssertIndexRange(face_orientation, 6);
729  (void)face;
730 
731  static const std::array<std::array<unsigned int, 3>, 6> table = {
732  {{{0, 2, 1}},
733  {{0, 1, 2}},
734  {{1, 2, 0}},
735  {{1, 0, 2}},
736  {{2, 1, 0}},
737  {{2, 0, 1}}}};
738 
739  return table[face_orientation][vertex];
740  }
741 
743  face_reference_cell_type(const unsigned int face_no) const override
744  {
745  (void)face_no;
746 
747  AssertIndexRange(face_no, n_faces());
748 
750  }
751 
752  unsigned int
753  face_to_cell_lines(const unsigned int face,
754  const unsigned int line,
755  const unsigned char face_orientation) const override
756  {
757  AssertIndexRange(face, n_faces());
758 
759  const static std::array<std::array<unsigned int, 3>, 4> table = {
760  {{{0, 1, 2}}, {{0, 3, 4}}, {{2, 5, 3}}, {{1, 4, 5}}}};
761 
762  return table[face][standard_to_real_face_line(
763  line, face, face_orientation)];
764  }
765 
766  unsigned int
768  const unsigned int face,
769  const unsigned int vertex,
770  const unsigned char face_orientation) const override
771  {
772  static const std::array<std::array<unsigned int, 3>, 4> table = {
773  {{{0, 1, 2}}, {{1, 0, 3}}, {{0, 2, 3}}, {{2, 1, 3}}}};
774 
775  return table[face][standard_to_real_face_vertex(
776  vertex, face, face_orientation)];
777  }
778 
779  virtual unsigned int
781  const unsigned int vertex_n) const override
782  {
783  AssertIndexRange(vertex_n, n_vertices());
784  return vertex_n;
785  }
786 
787  virtual unsigned int
788  exodusii_face_to_deal_face(const unsigned int face_n) const override
789  {
790  AssertIndexRange(face_n, n_faces());
791  constexpr std::array<unsigned int, 4> exodus_to_deal{{1, 3, 2, 0}};
792  return exodus_to_deal[face_n];
793  }
794  };
795 
796 
797 
801  struct Pyramid : public TensorProductBase<3>
802  {
803  unsigned int
804  n_vertices() const override
805  {
806  return 5;
807  }
808 
809  unsigned int
810  n_lines() const override
811  {
812  return 8;
813  }
814 
815  unsigned int
816  n_faces() const override
817  {
818  return 5;
819  }
820 
821  std::array<unsigned int, 2>
823  const unsigned int line) const override
824  {
825  Assert(false, ExcNotImplemented());
826 
827  static const std::array<unsigned int, 2> table[6] = {
828  {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}};
829 
830  return table[line];
831  }
832 
833  unsigned int
835  const unsigned int line,
836  const unsigned int face,
837  const unsigned char face_orientation) const override
838  {
839  Assert(false, ExcNotImplemented());
840 
841  (void)face;
842 
843  static const std::array<std::array<unsigned int, 3>, 6> table = {
844  {{{2, 1, 0}},
845  {{0, 1, 2}},
846  {{1, 2, 0}},
847  {{0, 2, 1}},
848  {{1, 0, 2}},
849  {{2, 0, 1}}}};
850 
851  return table[face_orientation][line];
852  }
853 
854  bool
856  const unsigned int line,
857  const unsigned char face_orientation_raw,
858  const unsigned char line_orientation) const override
859  {
860  (void)line;
861  (void)face_orientation_raw;
862 
863  return line_orientation;
864  }
865 
866  std::array<unsigned int, 2>
868  const unsigned int vertex) const override
869  {
870  static const std::array<unsigned int, 2> table[5] = {
871  {{0, 0}}, {{0, 1}}, {{0, 2}}, {{0, 3}}, {{1, 2}}};
872 
873  return table[vertex];
874  }
875 
876  unsigned int
878  const unsigned int vertex,
879  const unsigned int face,
880  const unsigned char face_orientation) const override
881  {
882  if (face == 0) // Quad
883  {
885  vertex,
886  get_bit(face_orientation, 0),
887  get_bit(face_orientation, 2),
888  get_bit(face_orientation, 1));
889  }
890  else // Tri
891  {
892  static const std::array<std::array<unsigned int, 3>, 6> table = {
893  {{{0, 2, 1}},
894  {{0, 1, 2}},
895  {{1, 2, 0}},
896  {{1, 0, 2}},
897  {{2, 1, 0}},
898  {{2, 0, 1}}}};
899 
900  return table[face_orientation][vertex];
901  }
902  }
903 
905  face_reference_cell_type(const unsigned int face_no) const override
906  {
907  AssertIndexRange(face_no, n_faces());
908 
909  if (face_no == 0)
911  else
913  }
914 
915  unsigned int
917  const unsigned int face,
918  const unsigned int vertex,
919  const unsigned char face_orientation) const override
920  {
921  AssertIndexRange(face, n_faces());
922  if (face == 0)
923  {
924  AssertIndexRange(vertex, 4);
925  }
926  else
927  {
928  AssertIndexRange(vertex, 3);
929  }
930  constexpr auto X = numbers::invalid_unsigned_int;
931  static const std::array<std::array<unsigned int, 4>, 5> table = {
932  {{{0, 1, 2, 3}},
933  {{0, 2, 4, X}},
934  {{3, 1, 4, X}},
935  {{1, 0, 4, X}},
936  {{2, 3, 4, X}}}};
937 
938  return table[face][standard_to_real_face_vertex(
939  vertex, face, face_orientation)];
940  }
941 
942  virtual unsigned int
944  const unsigned int vertex_n) const override
945  {
946  AssertIndexRange(vertex_n, n_vertices());
947  constexpr std::array<unsigned int, 5> exodus_to_deal{{0, 1, 3, 2, 4}};
948  return exodus_to_deal[vertex_n];
949  }
950 
951  virtual unsigned int
952  exodusii_face_to_deal_face(const unsigned int face_n) const override
953  {
954  AssertIndexRange(face_n, n_faces());
955  constexpr std::array<unsigned int, 5> exodus_to_deal{{3, 2, 4, 1, 0}};
956  return exodus_to_deal[face_n];
957  }
958  };
959 
960 
961 
965  struct Wedge : public TensorProductBase<3>
966  {
967  unsigned int
968  n_vertices() const override
969  {
970  return 6;
971  }
972 
973  unsigned int
974  n_lines() const override
975  {
976  return 9;
977  }
978 
979  unsigned int
980  n_faces() const override
981  {
982  return 5;
983  }
984 
985  std::array<unsigned int, 2>
987  const unsigned int line) const override
988  {
989  Assert(false, ExcNotImplemented());
990 
991  static const std::array<unsigned int, 2> table[6] = {
992  {{0, 0}}, {{0, 1}}, {{0, 2}}, {{1, 1}}, {{1, 2}}, {{2, 1}}};
993 
994  return table[line];
995  }
996 
997  unsigned int
999  const unsigned int line,
1000  const unsigned int face,
1001  const unsigned char face_orientation) const override
1002  {
1003  Assert(false, ExcNotImplemented());
1004 
1005  (void)face;
1006 
1007  static const std::array<std::array<unsigned int, 3>, 6> table = {
1008  {{{2, 1, 0}},
1009  {{0, 1, 2}},
1010  {{1, 2, 0}},
1011  {{0, 2, 1}},
1012  {{1, 0, 2}},
1013  {{2, 0, 1}}}};
1014 
1015  return table[face_orientation][line];
1016  }
1017 
1018  bool
1020  const unsigned int line,
1021  const unsigned char face_orientation_raw,
1022  const unsigned char line_orientation) const override
1023  {
1024  (void)line;
1025  (void)face_orientation_raw;
1026 
1027  return line_orientation;
1028  }
1029 
1030  std::array<unsigned int, 2>
1032  const unsigned int vertex) const override
1033  {
1034  static const std::array<std::array<unsigned int, 2>, 6> table = {
1035  {{{0, 1}}, {{0, 0}}, {{0, 2}}, {{1, 0}}, {{1, 1}}, {{1, 2}}}};
1036 
1037  return table[vertex];
1038  }
1039 
1040  unsigned int
1042  const unsigned int vertex,
1043  const unsigned int face,
1044  const unsigned char face_orientation) const override
1045  {
1046  if (face > 1) // QUAD
1047  {
1049  vertex,
1050  get_bit(face_orientation, 0),
1051  get_bit(face_orientation, 2),
1052  get_bit(face_orientation, 1));
1053  }
1054  else // TRI
1055  {
1056  static const std::array<std::array<unsigned int, 3>, 6> table = {
1057  {{{0, 2, 1}},
1058  {{0, 1, 2}},
1059  {{1, 2, 0}},
1060  {{1, 0, 2}},
1061  {{2, 1, 0}},
1062  {{2, 0, 1}}}};
1063 
1064  return table[face_orientation][vertex];
1065  }
1066  }
1067 
1069  face_reference_cell_type(const unsigned int face_no) const override
1070  {
1071  AssertIndexRange(face_no, n_faces());
1072 
1073  if (face_no > 1)
1075  else
1076  return ReferenceCell::Type::Tri;
1077  }
1078 
1079  unsigned int
1081  const unsigned int face,
1082  const unsigned int vertex,
1083  const unsigned char face_orientation) const override
1084  {
1085  AssertIndexRange(face, n_faces());
1086  if (face < 2)
1087  {
1088  AssertIndexRange(vertex, 3);
1089  }
1090  else
1091  {
1092  AssertIndexRange(vertex, 4);
1093  }
1094  constexpr auto X = numbers::invalid_unsigned_int;
1095  static const std::array<std::array<unsigned int, 4>, 6> table = {
1096  {{{1, 0, 2, X}},
1097  {{3, 4, 5, X}},
1098  {{0, 1, 3, 4}},
1099  {{1, 2, 4, 5}},
1100  {{2, 0, 5, 3}}}};
1101 
1102  return table[face][standard_to_real_face_vertex(
1103  vertex, face, face_orientation)];
1104  }
1105 
1106  virtual unsigned int
1108  const unsigned int vertex_n) const override
1109  {
1110  AssertIndexRange(vertex_n, n_vertices());
1111  constexpr std::array<unsigned int, 6> exodus_to_deal{
1112  {2, 1, 0, 5, 4, 3}};
1113  return exodus_to_deal[vertex_n];
1114  }
1115 
1116  virtual unsigned int
1117  exodusii_face_to_deal_face(const unsigned int face_n) const override
1118  {
1119  AssertIndexRange(face_n, n_faces());
1120  constexpr std::array<unsigned int, 6> exodus_to_deal{{3, 4, 2, 0, 1}};
1121  return exodus_to_deal[face_n];
1122  }
1123  };
1124 
1125 
1126 
1130  struct Hex : public TensorProductBase<3>
1131  {
1132  std::array<unsigned int, 2>
1134  const unsigned int line) const override
1135  {
1137  }
1138 
1139  unsigned int
1141  const unsigned int line,
1142  const unsigned int face,
1143  const unsigned char face_orientation) const override
1144  {
1145  (void)face;
1146 
1148  line,
1149  get_bit(face_orientation, 0),
1150  get_bit(face_orientation, 2),
1151  get_bit(face_orientation, 1));
1152  }
1153 
1154  bool
1156  const unsigned int line,
1157  const unsigned char face_orientation_raw,
1158  const unsigned char line_orientation) const override
1159  {
1160  static const bool bool_table[2][2][2][2] = {
1161  {{{true, false}, // lines 0/1, face_orientation=false,
1162  // face_flip=false, face_rotation=false and true
1163  {false, true}}, // lines 0/1, face_orientation=false,
1164  // face_flip=true, face_rotation=false and true
1165  {{true, true}, // lines 0/1, face_orientation=true,
1166  // face_flip=false, face_rotation=false and true
1167  {false, false}}}, // lines 0/1, face_orientation=true,
1168  // face_flip=true, face_rotation=false and true
1169 
1170  {{{true, true}, // lines 2/3 ...
1171  {false, false}},
1172  {{true, false}, {false, true}}}};
1173 
1174  const bool face_orientation = get_bit(face_orientation_raw, 0);
1175  const bool face_flip = get_bit(face_orientation_raw, 2);
1176  const bool face_rotation = get_bit(face_orientation_raw, 1);
1177 
1178  return (
1179  static_cast<bool>(line_orientation) ==
1180  bool_table[line / 2][face_orientation][face_flip][face_rotation]);
1181  }
1182 
1183  std::array<unsigned int, 2>
1185  const unsigned int vertex) const override
1186  {
1188  vertex);
1189  }
1190 
1191  unsigned int
1193  const unsigned int vertex,
1194  const unsigned int face,
1195  const unsigned char face_orientation) const override
1196  {
1197  (void)face;
1198 
1200  vertex,
1201  get_bit(face_orientation, 0),
1202  get_bit(face_orientation, 2),
1203  get_bit(face_orientation, 1));
1204  }
1205 
1207  face_reference_cell_type(const unsigned int face_no) const override
1208  {
1209  (void)face_no;
1211  }
1212 
1213  virtual unsigned int
1215  const unsigned int vertex_n) const override
1216  {
1217  AssertIndexRange(vertex_n, n_vertices());
1218  constexpr std::array<unsigned int, 8> exodus_to_deal{
1219  {0, 1, 3, 2, 4, 5, 7, 6}};
1220  return exodus_to_deal[vertex_n];
1221  }
1222 
1223  virtual unsigned int
1224  exodusii_face_to_deal_face(const unsigned int face_n) const override
1225  {
1226  AssertIndexRange(face_n, n_faces());
1227  constexpr std::array<unsigned int, 6> exodus_to_deal{
1228  {2, 1, 3, 0, 4, 5}};
1229  return exodus_to_deal[face_n];
1230  }
1231  };
1232 
1238  {
1239  static const std::
1240  array<std::unique_ptr<ReferenceCell::internal::Info::Base>, 8>
1241  gei{{std::make_unique<ReferenceCell::internal::Info::Vertex>(),
1242  std::make_unique<ReferenceCell::internal::Info::Line>(),
1243  std::make_unique<ReferenceCell::internal::Info::Tri>(),
1244  std::make_unique<ReferenceCell::internal::Info::Quad>(),
1245  std::make_unique<ReferenceCell::internal::Info::Tet>(),
1246  std::make_unique<ReferenceCell::internal::Info::Pyramid>(),
1247  std::make_unique<ReferenceCell::internal::Info::Wedge>(),
1248  std::make_unique<ReferenceCell::internal::Info::Hex>()}};
1249  AssertIndexRange(static_cast<std::uint8_t>(type), 8);
1250  return *gei[static_cast<std::uint8_t>(type)];
1251  }
1252 
1258  get_face(const ReferenceCell::Type &type, const unsigned int face_no)
1259  {
1260  return get_cell(get_cell(type).face_reference_cell_type(face_no));
1261  }
1262 
1263  } // namespace Info
1264  } // namespace internal
1265 } // namespace ReferenceCell
1266 
1267 
1269 
1270 #endif
unsigned int n_faces() const override
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const override
virtual unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const override
static bool get_bit(const unsigned char number, const unsigned int n)
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const override
static const unsigned int invalid_unsigned_int
Definition: types.h:196
virtual unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const override
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
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char line_orientation) const override
ReferenceCell::Type face_reference_cell_type(const unsigned int face_no) const override
virtual unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const override
std::array< unsigned int, 2 > standard_line_to_face_and_line_index(const unsigned int line) const override
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)
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const override
unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const unsigned char face_orientation) const override
unsigned int n_lines() const override
virtual unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) 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)
ReferenceCell::Type face_reference_cell_type(const unsigned int face_no) const override
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const override
virtual unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const override
Type get_hypercube(const unsigned int dim)
unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const unsigned char face_orientation) const override
virtual std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const
unsigned int standard_to_real_face_line(const unsigned int line, const unsigned int face, const unsigned char face_orientation) const override
ReferenceCell::Type face_reference_cell_type(const unsigned int face_no) const override
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const override
const ReferenceCell::internal::Info::Base & get_cell(const ReferenceCell::Type &type)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1691
virtual unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const override
virtual unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const override
ReferenceCell::Type face_reference_cell_type(const unsigned int face_no) const override
virtual unsigned int n_faces() const
virtual unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const override
unsigned int n_faces() const override
std::array< unsigned int, 2 > standard_line_to_face_and_line_index(const unsigned int line) const override
static const char U
virtual unsigned int standard_to_real_face_line(const unsigned int line, const unsigned int face, const unsigned char face_orientation) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
unsigned int n_vertices() const override
ReferenceCell::Type face_reference_cell_type(const unsigned int face_no) const override
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const override
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
bool combine_face_and_line_orientation(const unsigned int line, const unsigned char face_orientation_raw, const unsigned char line_orientation) const override
std::array< unsigned int, 2 > standard_line_to_face_and_line_index(const unsigned int line) const override
virtual unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char face_orientation) const
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char face_orientation) const override
static std::array< unsigned int, 2 > standard_quad_vertex_to_line_vertex_index(const unsigned int vertex)
virtual unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const override
const ReferenceCell::internal::Info::Base & get_face(const ReferenceCell::Type &type, const unsigned int face_no)
ReferenceCell::Type face_reference_cell_type(const unsigned int face_no) const override
ReferenceCell::Type face_reference_cell_type(const unsigned int face_no) const override
static ::ExceptionBase & ExcMessage(std::string arg1)
virtual ReferenceCell::Type face_reference_cell_type(const unsigned int face_no) const
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const override
virtual unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const
#define Assert(cond, exc)
Definition: exceptions.h:1466
unsigned int n_lines() const override
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char face_orientation) const override
virtual bool combine_face_and_line_orientation(const unsigned int line, const unsigned char face_orientation, const unsigned char line_orientation) const
static void set_bit(unsigned char &number, const unsigned int n, const bool x)
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const override
bool combine_face_and_line_orientation(const unsigned int line, const unsigned char face_orientation_raw, const unsigned char line_orientation) const override
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const override
virtual unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const override
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:372
static unsigned int standard_to_real_line_vertex(const unsigned int vertex, const bool line_orientation=true)
virtual unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const override
std::string to_string(const T &t)
Definition: patterns.h:2342
virtual unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const override
unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const unsigned char face_orientation) const override
unsigned int n_faces() const override
Type get_simplex(const unsigned int dim)
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)
virtual unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const unsigned char face_orientation) const
virtual unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const override
virtual std::array< unsigned int, 2 > standard_line_to_face_and_line_index(const unsigned int line) const
unsigned int n_vertices() const override
unsigned int n_lines() const override
unsigned int n_lines() const override
virtual unsigned int n_lines() const
virtual unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const override
unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const unsigned char face_orientation) const override
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:371
virtual unsigned int n_vertices() const
virtual unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const
std::array< unsigned int, 2 > standard_vertex_to_face_and_vertex_index(const unsigned int vertex) const override
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char face_orientation) const override
Type n_vertices_to_type(const int dim, const unsigned int n_vertices)
bool combine_face_and_line_orientation(const unsigned int line, const unsigned char face_orientation_raw, const unsigned char line_orientation) const override
static std::array< unsigned int, 2 > standard_hex_vertex_to_quad_vertex_index(const unsigned int vertex)
static ::ExceptionBase & ExcNotImplemented()
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char face_orientation) const override
unsigned int n_vertices() const override
unsigned int n_vertices() const override
unsigned int standard_to_real_face_line(const unsigned int line, const unsigned int face, const unsigned char face_orientation) const override
std::array< unsigned int, 2 > standard_line_to_face_and_line_index(const unsigned int line) const override
unsigned int standard_to_real_face_line(const unsigned int line, const unsigned int face, const unsigned char face_orientation) const override
virtual unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const override
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
ReferenceCell::Type face_reference_cell_type(const unsigned int face_no) const override
unsigned int standard_to_real_face_line(const unsigned int line, const unsigned int face, const unsigned char face_orientation) const override
bool combine_face_and_line_orientation(const unsigned int line, const unsigned char face_orientation_raw, const unsigned char line_orientation) const override
static std::array< unsigned int, 2 > standard_hex_line_to_quad_line_index(const unsigned int line)
unsigned int standard_to_real_face_vertex(const unsigned int vertex, const unsigned int face, const unsigned char line_orientation) const override
unsigned int n_faces() const override
virtual unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const override