Reference documentation for deal.II version GIT 29f9da0a34 2023-12-07 10:00: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\}}\)
geometry_info.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2023 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_geometry_info_h
17 #define dealii_geometry_info_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 #include <deal.II/base/ndarray.h>
24 #include <deal.II/base/point.h>
26 
27 #include <array>
28 #include <cstdint>
29 
30 
31 
33 
34 #ifndef DOXYGEN
35 namespace internal
36 {
37  namespace GeometryInfoHelper
38  {
39  // A struct that holds the values for all the arrays we want to initialize
40  // in GeometryInfo
41  template <int dim>
42  struct Initializers;
43 
44  template <>
45  struct Initializers<1>
46  {
47  static constexpr std::array<unsigned int, 2>
48  ucd_to_deal()
49  {
50  return {{0, 1}};
51  }
52 
53  static constexpr std::array<unsigned int, 2>
54  unit_normal_direction()
55  {
56  return {{0, 0}};
57  }
58 
59  static constexpr std::array<int, 2>
60  unit_normal_orientation()
61  {
62  return {{-1, 1}};
63  }
64 
65  static constexpr std::array<Tensor<1, 1>, 2>
66  unit_normal_vector()
67  {
68  return {{Tensor<1, 1>{{-1}}, Tensor<1, 1>{{1}}}};
69  }
70 
71  static constexpr ::ndarray<Tensor<1, 1>, 2, 0>
72  unit_tangential_vectors()
73  {
74  return {{{{}}, {{}}}};
75  }
76 
77  static constexpr std::array<unsigned int, 2>
78  opposite_face()
79  {
80  return {{1, 0}};
81  }
82 
83  static constexpr std::array<unsigned int, 2>
84  dx_to_deal()
85  {
86  return {{0, 1}};
87  }
88 
89  static constexpr ::ndarray<unsigned int, 2, 1>
90  vertex_to_face()
91  {
92  return {{{{0}}, {{1}}}};
93  }
94  };
95 
96  template <>
97  struct Initializers<2>
98  {
99  static constexpr std::array<unsigned int, 4>
100  ucd_to_deal()
101  {
102  return {{0, 1, 3, 2}};
103  }
104 
105  static constexpr std::array<unsigned int, 4>
106  unit_normal_direction()
107  {
108  return {{0, 0, 1, 1}};
109  }
110 
111  static constexpr std::array<int, 4>
112  unit_normal_orientation()
113  {
114  return {{-1, 1, -1, 1}};
115  }
116 
117  static constexpr std::array<Tensor<1, 2>, 4>
118  unit_normal_vector()
119  {
120  return {{Tensor<1, 2>{{-1., 0.}},
121  Tensor<1, 2>{{1., 0.}},
122  Tensor<1, 2>{{0., -1.}},
123  Tensor<1, 2>{{0., 1.}}}};
124  }
125 
126  static constexpr ::ndarray<Tensor<1, 2>, 4, 1>
127  unit_tangential_vectors()
128  {
129  return {{{{Tensor<1, 2>{{0, -1}}}},
130  {{Tensor<1, 2>{{0, 1}}}},
131  {{Tensor<1, 2>{{1, 0}}}},
132  {{Tensor<1, 2>{{-1, 0}}}}}};
133  }
134 
135  static constexpr std::array<unsigned int, 4>
136  opposite_face()
137  {
138  return {{1, 0, 3, 2}};
139  }
140 
141  static constexpr std::array<unsigned int, 4>
142  dx_to_deal()
143  {
144  return {{0, 2, 1, 3}};
145  }
146 
147  static constexpr ::ndarray<unsigned int, 4, 2>
148  vertex_to_face()
149  {
150  return {{{{0, 2}}, {{1, 2}}, {{0, 3}}, {{1, 3}}}};
151  }
152  };
153 
154  template <>
155  struct Initializers<3>
156  {
157  static constexpr std::array<unsigned int, 8>
158  ucd_to_deal()
159  {
160  return {{0, 4, 5, 1, 2, 6, 7, 3}};
161  }
162 
163  static constexpr std::array<unsigned int, 6>
164  unit_normal_direction()
165  {
166  return {{0, 0, 1, 1, 2, 2}};
167  }
168 
169  static constexpr std::array<int, 6>
170  unit_normal_orientation()
171  {
172  return {{-1, 1, -1, 1, -1, 1}};
173  }
174 
175  static constexpr std::array<Tensor<1, 3>, 6>
176  unit_normal_vector()
177  {
178  return {{Tensor<1, 3>{{-1, 0, 0}},
179  Tensor<1, 3>{{1, 0, 0}},
180  Tensor<1, 3>{{0, -1, 0}},
181  Tensor<1, 3>{{0, 1, 0}},
182  Tensor<1, 3>{{0, 0, -1}},
183  Tensor<1, 3>{{0, 0, 1}}}};
184  }
185 
186  static constexpr ::ndarray<Tensor<1, 3>, 6, 2>
187  unit_tangential_vectors()
188  {
189  return {{{{Tensor<1, 3>{{0, -1, 0}}, Tensor<1, 3>{{0, 0, 1}}}},
190  {{Tensor<1, 3>{{0, 1, 0}}, Tensor<1, 3>{{0, 0, 1}}}},
191  {{Tensor<1, 3>{{0, 0, -1}}, Tensor<1, 3>{{1, 0, 0}}}},
192  {{Tensor<1, 3>{{0, 0, 1}}, Tensor<1, 3>{{1, 0, 0}}}},
193  {{Tensor<1, 3>{{-1, 0, 0}}, Tensor<1, 3>{{0, 1, 0}}}},
194  {{Tensor<1, 3>{{1, 0, 0}}, Tensor<1, 3>{{0, 1, 0}}}}}};
195  }
196 
197  static constexpr std::array<unsigned int, 6>
198  opposite_face()
199  {
200  return {{1, 0, 3, 2, 5, 4}};
201  }
202 
203  static constexpr std::array<unsigned int, 8>
204  dx_to_deal()
205  {
206  return {{0, 4, 2, 6, 1, 5, 3, 7}};
207  }
208 
209  static constexpr ::ndarray<unsigned int, 8, 3>
210  vertex_to_face()
211  {
212  return {{{{0, 2, 4}},
213  {{1, 2, 4}},
214  {{0, 3, 4}},
215  {{1, 3, 4}},
216  {{0, 2, 5}},
217  {{1, 2, 5}},
218  {{0, 3, 5}},
219  {{1, 3, 5}}}};
220  }
221  };
222 
223  template <>
224  struct Initializers<4>
225  {
226  static constexpr std::array<unsigned int, 16>
227  ucd_to_deal()
228  {
245  }
246 
247  static constexpr std::array<unsigned int, 8>
248  unit_normal_direction()
249  {
250  return {{0, 0, 1, 1, 2, 2, 3, 3}};
251  }
252 
253  static constexpr std::array<int, 8>
254  unit_normal_orientation()
255  {
256  return {{-1, 1, -1, 1, -1, 1, -1, 1}};
257  }
258 
259  static constexpr std::array<Tensor<1, 4>, 8>
260  unit_normal_vector()
261  {
262  return {{Tensor<1, 4>{{-1, 0, 0, 0}},
263  Tensor<1, 4>{{1, 0, 0, 0}},
264  Tensor<1, 4>{{0, -1, 0, 0}},
265  Tensor<1, 4>{{0, 1, 0, 0}},
266  Tensor<1, 4>{{0, 0, -1, 0}},
267  Tensor<1, 4>{{0, 0, 1, 0}},
268  Tensor<1, 4>{{0, 0, 0, -1}},
269  Tensor<1, 4>{{0, 0, 0, 1}}}};
270  }
271 
272  static constexpr ::ndarray<Tensor<1, 4>, 8, 3>
273  unit_tangential_vectors()
274  {
275  return {{{{Tensor<1, 4>{{0, -1, 0, 0}},
276  Tensor<1, 4>{{0, 0, 1, 0}},
277  Tensor<1, 4>{{0, 0, 0, 1}}}},
278  {{Tensor<1, 4>{{0, 1, 0, 0}},
279  Tensor<1, 4>{{0, 0, 1, 0}},
280  Tensor<1, 4>{{0, 0, 0, 1}}}},
281  {{Tensor<1, 4>{{0, 0, -1, 0}},
282  Tensor<1, 4>{{0, 0, 0, 1}},
283  Tensor<1, 4>{{1, 0, 0, 0}}}},
284  {{Tensor<1, 4>{{0, 0, 1, 0}},
285  Tensor<1, 4>{{0, 0, 0, 1}},
286  Tensor<1, 4>{{1, 0, 0, 0}}}},
287  {{Tensor<1, 4>{{0, 0, 0, -1}},
288  Tensor<1, 4>{{1, 0, 0, 0}},
289  Tensor<1, 4>{{0, 1, 0, 0}}}},
290  {{Tensor<1, 4>{{0, 0, 0, 1}},
291  Tensor<1, 4>{{1, 0, 0, 0}},
292  Tensor<1, 4>{{0, 1, 0, 0}}}},
293  {{Tensor<1, 4>{{-1, 0, 0, 0}},
294  Tensor<1, 4>{{0, 1, 0, 0}},
295  Tensor<1, 4>{{0, 0, 1, 0}}}},
296  {{Tensor<1, 4>{{1, 0, 0, 0}},
297  Tensor<1, 4>{{0, 1, 0, 0}},
298  Tensor<1, 4>{{0, 0, 1, 0}}}}}};
299  }
300 
301  static constexpr std::array<unsigned int, 8>
302  opposite_face()
303  {
304  return {{1, 0, 3, 2, 5, 4, 7, 6}};
305  }
306 
307  static constexpr std::array<unsigned int, 16>
308  dx_to_deal()
309  {
326  }
327 
328  static constexpr ::ndarray<unsigned int, 16, 4>
329  vertex_to_face()
330  {
395  }
396  };
397  } // namespace GeometryInfoHelper
398 } // namespace internal
399 #endif // DOXYGEN
400 
401 
418 {
419 public:
426  enum Object
427  {
431  vertex = 0,
435  line = 1,
439  quad = 2,
443  hex = 3
444  };
445 
450  GeometryPrimitive(const Object object);
451 
457  GeometryPrimitive(const unsigned int object_dimension);
458 
463  operator unsigned int() const;
464 
465 private:
470 };
471 
472 
473 
486 template <int dim>
488 {
523  enum Possibilities : std::uint8_t
524  {
529 
541  isotropic_refinement = 0xFF
542  };
543 };
544 
545 
546 
556 template <>
558 {
593  enum Possibilities : std::uint8_t
594  {
602  cut_x = 1,
606  isotropic_refinement = cut_x
607  };
608 };
609 
610 
611 
622 template <>
624 {
659  enum Possibilities : std::uint8_t
660  {
668  cut_x = 1,
672  cut_y = 2,
676  cut_xy = cut_x | cut_y,
677 
681  isotropic_refinement = cut_xy
682  };
683 };
684 
685 
686 
697 template <>
699 {
734  enum Possibilities : std::uint8_t
735  {
743  cut_x = 1,
747  cut_y = 2,
751  cut_xy = cut_x | cut_y,
755  cut_z = 4,
759  cut_xz = cut_x | cut_z,
763  cut_yz = cut_y | cut_z,
767  cut_xyz = cut_x | cut_y | cut_z,
768 
772  isotropic_refinement = cut_xyz
773  };
774 };
775 
776 
777 
788 template <int dim>
790 {
791 public:
795  static constexpr unsigned int n_refinement_cases = (1 << dim);
796 
801 
807  const typename RefinementPossibilities<dim>::Possibilities refinement_case);
808 
814  explicit RefinementCase(const std::uint8_t refinement_case);
815 
827  DEAL_II_HOST_DEVICE operator std::uint8_t() const;
828 
834  operator|(const RefinementCase &r) const;
835 
841  operator&(const RefinementCase &r) const;
842 
851  operator~() const;
852 
853 
859  static RefinementCase
860  cut_axis(const unsigned int i);
861 
870  static std::array<RefinementCase<dim>, n_refinement_cases>
872 
876  static std::size_t
878 
884  template <class Archive>
885  void
886  serialize(Archive &ar, const unsigned int version);
887 
893  int,
894  << "The refinement flags given (" << arg1
895  << ") contain set bits that do not "
896  << "make sense for the space dimension of the object to which they are applied.");
897 
898 private:
903  std::uint8_t value : (dim > 0 ? dim : 1);
904 };
905 
906 
907 
908 namespace internal
909 {
929  template <int dim>
931  {
936  {
941 
945  case_isotropic = static_cast<std::uint8_t>(-1)
946  };
947  };
948 
949 
958  template <>
960  {
967  {
972 
977  };
978  };
979 
980 
981 
991  template <>
993  {
1000  {
1005 
1010  };
1011  };
1012 
1013 
1014 
1025  template <>
1027  {
1035  {
1043  case_x = 1,
1047  case_isotropic = case_x
1048  };
1049  };
1050 
1051 
1052 
1144  template <>
1146  {
1154  {
1156  case_x = 1,
1157  case_x1y = 2,
1158  case_x2y = 3,
1159  case_x1y2y = 4,
1160  case_y = 5,
1161  case_y1x = 6,
1162  case_y2x = 7,
1163  case_y1x2x = 8,
1164  case_xy = 9,
1165 
1166  case_isotropic = case_xy
1167  };
1168  };
1169 
1170 
1171 
1178  template <int dim>
1179  class SubfaceCase : public SubfacePossibilities<dim>
1180  {
1181  public:
1188  subface_possibility);
1189 
1201  operator std::uint8_t() const;
1202 
1206  static constexpr std::size_t
1208 
1214  int,
1215  << "The subface case given (" << arg1 << ") does not make sense "
1216  << "for the space dimension of the object to which they are applied.");
1217 
1218  private:
1223  std::uint8_t value : (dim == 3 ? 4 : 1);
1224  };
1225 
1226 } // namespace internal
1227 
1228 
1229 
1230 template <int dim>
1231 struct GeometryInfo;
1232 
1233 
1234 
1254 template <>
1255 struct GeometryInfo<0>
1256 {
1264  static constexpr unsigned int max_children_per_cell = 1;
1265 
1269  static constexpr unsigned int faces_per_cell = 0;
1270 
1287  static std::array<unsigned int, 0>
1289 
1297  static constexpr unsigned int max_children_per_face = 0;
1298 
1304  static unsigned int
1305  n_children(const RefinementCase<0> &refinement_case);
1306 
1310  static constexpr unsigned int vertices_per_cell = 1;
1311 
1330  static std::array<unsigned int, vertices_per_cell>
1332 
1356  static unsigned int
1357  face_to_cell_vertices(const unsigned int face,
1358  const unsigned int vertex,
1359  const bool face_orientation = true,
1360  const bool face_flip = false,
1361  const bool face_rotation = false);
1362 
1377  static unsigned int
1378  face_to_cell_lines(const unsigned int face,
1379  const unsigned int line,
1380  const bool face_orientation = true,
1381  const bool face_flip = false,
1382  const bool face_rotation = false);
1383 
1390  static constexpr unsigned int vertices_per_face = 0;
1391 
1395  static constexpr unsigned int lines_per_face = 0;
1396 
1400  static constexpr unsigned int quads_per_face = 0;
1401 
1405  static constexpr unsigned int lines_per_cell = 0;
1406 
1410  static constexpr unsigned int quads_per_cell = 0;
1411 
1415  static constexpr unsigned int hexes_per_cell = 0;
1416 
1434  static const std::array<unsigned int, vertices_per_cell> ucd_to_deal;
1435 
1449  static const std::array<unsigned int, vertices_per_cell> dx_to_deal;
1450 };
1451 
1452 
1453 
1985 template <int dim>
1987 {
1995  static constexpr unsigned int max_children_per_cell = 1 << dim;
1996 
2000  static constexpr unsigned int faces_per_cell = 2 * dim;
2001 
2020 
2028  static constexpr unsigned int max_children_per_face =
2030 
2034  static constexpr unsigned int vertices_per_cell = 1 << dim;
2035 
2053 
2057  static constexpr unsigned int vertices_per_face =
2059 
2063  static constexpr unsigned int lines_per_face =
2065 
2069  static constexpr unsigned int quads_per_face =
2071 
2081  static constexpr unsigned int lines_per_cell =
2082  (2 * GeometryInfo<dim - 1>::lines_per_cell +
2084 
2092  static constexpr unsigned int quads_per_cell =
2093  (2 * GeometryInfo<dim - 1>::quads_per_cell +
2094  GeometryInfo<dim - 1>::lines_per_cell);
2095 
2099  static constexpr unsigned int hexes_per_cell =
2100  (2 * GeometryInfo<dim - 1>::hexes_per_cell +
2101  GeometryInfo<dim - 1>::quads_per_cell);
2102 
2120  static constexpr std::array<unsigned int, vertices_per_cell> ucd_to_deal =
2121  internal::GeometryInfoHelper::Initializers<dim>::ucd_to_deal();
2122 
2136  static constexpr std::array<unsigned int, vertices_per_cell> dx_to_deal =
2137  internal::GeometryInfoHelper::Initializers<dim>::dx_to_deal();
2138 
2151  internal::GeometryInfoHelper::Initializers<dim>::vertex_to_face();
2152 
2157  static unsigned int
2158  n_children(const RefinementCase<dim> &refinement_case);
2159 
2164  static unsigned int
2166 
2176  static double
2178  const unsigned int subface_no);
2179 
2185  static RefinementCase<dim - 1>
2186  face_refinement_case(const RefinementCase<dim> &cell_refinement_case,
2187  const unsigned int face_no,
2188  const bool face_orientation = true,
2189  const bool face_flip = false,
2190  const bool face_rotation = false);
2191 
2197  static RefinementCase<dim>
2200  const unsigned int face_no,
2201  const bool face_orientation = true,
2202  const bool face_flip = false,
2203  const bool face_rotation = false);
2204 
2209  static RefinementCase<1>
2210  line_refinement_case(const RefinementCase<dim> &cell_refinement_case,
2211  const unsigned int line_no);
2212 
2217  static RefinementCase<dim>
2219 
2266  static unsigned int
2268  const unsigned int face,
2269  const unsigned int subface,
2270  const bool face_orientation = true,
2271  const bool face_flip = false,
2272  const bool face_rotation = false,
2275 
2289  static unsigned int
2290  line_to_cell_vertices(const unsigned int line, const unsigned int vertex);
2291 
2312  static unsigned int
2313  face_to_cell_vertices(const unsigned int face,
2314  const unsigned int vertex,
2315  const bool face_orientation = true,
2316  const bool face_flip = false,
2317  const bool face_rotation = false);
2318 
2330  static unsigned int
2331  face_to_cell_lines(const unsigned int face,
2332  const unsigned int line,
2333  const bool face_orientation = true,
2334  const bool face_flip = false,
2335  const bool face_rotation = false);
2336 
2346  static unsigned int
2347  standard_to_real_face_vertex(const unsigned int vertex,
2348  const bool face_orientation = true,
2349  const bool face_flip = false,
2350  const bool face_rotation = false);
2351 
2361  static unsigned int
2362  real_to_standard_face_vertex(const unsigned int vertex,
2363  const bool face_orientation = true,
2364  const bool face_flip = false,
2365  const bool face_rotation = false);
2366 
2376  static unsigned int
2377  standard_to_real_face_line(const unsigned int line,
2378  const bool face_orientation = true,
2379  const bool face_flip = false,
2380  const bool face_rotation = false);
2381 
2387  static unsigned int
2388  standard_to_real_line_vertex(const unsigned int vertex,
2389  const bool line_orientation = true);
2390 
2398  static std::array<unsigned int, 2>
2399  standard_quad_vertex_to_line_vertex_index(const unsigned int vertex);
2400 
2408  static std::array<unsigned int, 2>
2409  standard_hex_vertex_to_quad_vertex_index(const unsigned int vertex);
2410 
2418  static std::array<unsigned int, 2>
2419  standard_hex_line_to_quad_line_index(const unsigned int line);
2420 
2430  static unsigned int
2431  real_to_standard_face_line(const unsigned int line,
2432  const bool face_orientation = true,
2433  const bool face_flip = false,
2434  const bool face_rotation = false);
2435 
2441  static Point<dim>
2442  unit_cell_vertex(const unsigned int vertex);
2443 
2453  static unsigned int
2455 
2463  static Point<dim>
2465  const unsigned int child_index,
2466  const RefinementCase<dim> refine_case =
2468 
2474  static Point<dim>
2476  const unsigned int child_index,
2477  const RefinementCase<dim> refine_case =
2479 
2484  static bool
2486 
2498  static bool
2499  is_inside_unit_cell(const Point<dim> &p, const double eps);
2500 
2505  template <typename Number = double>
2506  static Point<dim, Number>
2508 
2514  static double
2516 
2521  static double
2522  d_linear_shape_function(const Point<dim> &xi, const unsigned int i);
2523 
2528  static Tensor<1, dim>
2529  d_linear_shape_function_gradient(const Point<dim> &xi, const unsigned int i);
2530 
2582  template <int spacedim>
2583  static void
2585 #ifndef DEAL_II_CXX14_CONSTEXPR_BUG
2587  Tensor<spacedim - dim, spacedim> (&forms)[vertices_per_cell]);
2588 #else
2589  (const Point<spacedim> *vertices, Tensor<spacedim - dim, spacedim> *forms);
2590 #endif
2591 
2601  static constexpr std::array<unsigned int, faces_per_cell>
2603  internal::GeometryInfoHelper::Initializers<dim>::unit_normal_direction();
2604 
2621  static constexpr std::array<int, faces_per_cell> unit_normal_orientation =
2622  internal::GeometryInfoHelper::Initializers<dim>::unit_normal_orientation();
2623 
2634  static constexpr std::array<Tensor<1, dim>, faces_per_cell>
2636  internal::GeometryInfoHelper::Initializers<dim>::unit_normal_vector();
2637 
2651  static constexpr ndarray<Tensor<1, dim>, faces_per_cell, dim - 1>
2652  unit_tangential_vectors = internal::GeometryInfoHelper::Initializers<
2654 
2660  static constexpr std::array<unsigned int, faces_per_cell> opposite_face =
2661  internal::GeometryInfoHelper::Initializers<dim>::opposite_face();
2662 
2663 
2668  double,
2669  << "The coordinates must satisfy 0 <= x_i <= 1, "
2670  << "but here we have x_i=" << arg1);
2671 
2676  int,
2677  int,
2678  int,
2679  << "RefinementCase<dim> " << arg1 << ": face " << arg2
2680  << " has no subface " << arg3);
2681 };
2682 
2683 
2684 
2685 #ifndef DOXYGEN
2686 
2687 
2688 /* -------------- declaration of explicit specializations ------------- */
2689 
2690 
2691 template <>
2694  const unsigned int i);
2695 template <>
2698  const unsigned int i);
2699 template <>
2702  const unsigned int i);
2703 
2704 
2705 
2706 /* -------------- inline functions ------------- */
2707 
2708 
2709 inline GeometryPrimitive::GeometryPrimitive(const Object object)
2710  : object(object)
2711 {}
2712 
2713 
2714 
2715 inline GeometryPrimitive::GeometryPrimitive(const unsigned int object_dimension)
2716  : object(static_cast<Object>(object_dimension))
2717 {}
2718 
2719 
2720 inline GeometryPrimitive::operator unsigned int() const
2721 {
2722  return static_cast<unsigned int>(object);
2723 }
2724 
2725 
2726 
2727 namespace internal
2728 {
2729  template <int dim>
2731  const typename SubfacePossibilities<dim>::Possibilities subface_possibility)
2732  : value(subface_possibility)
2733  {}
2734 
2735 
2736  template <int dim>
2737  inline SubfaceCase<dim>::operator std::uint8_t() const
2738  {
2739  return value;
2740  }
2741 
2742 
2743 } // namespace internal
2744 
2745 
2746 template <int dim>
2747 inline RefinementCase<dim>
2748 RefinementCase<dim>::cut_axis(const unsigned int)
2749 {
2750  Assert(false, ExcInternalError());
2751  return static_cast<std::uint8_t>(-1);
2752 }
2753 
2754 
2755 template <>
2756 inline RefinementCase<1>
2757 RefinementCase<1>::cut_axis(const unsigned int i)
2758 {
2759  AssertIndexRange(i, 1);
2760 
2762  return options[i];
2763 }
2764 
2765 
2766 
2767 template <>
2768 inline RefinementCase<2>
2769 RefinementCase<2>::cut_axis(const unsigned int i)
2770 {
2771  AssertIndexRange(i, 2);
2772 
2775  return options[i];
2776 }
2777 
2778 
2779 
2780 template <>
2781 inline RefinementCase<3>
2782 RefinementCase<3>::cut_axis(const unsigned int i)
2783 {
2784  AssertIndexRange(i, 3);
2785 
2789  return options[i];
2790 }
2791 
2792 
2793 
2794 template <>
2795 inline std::array<RefinementCase<1>, 2>
2797 {
2800 }
2801 
2802 
2803 
2804 template <>
2805 inline std::array<RefinementCase<2>, 4>
2807 {
2812 }
2813 
2814 
2815 
2816 template <>
2817 inline std::array<RefinementCase<3>, 8>
2819 {
2828 }
2829 
2830 
2831 
2832 template <int dim>
2834  : value(RefinementPossibilities<dim>::no_refinement)
2835 {}
2836 
2837 
2838 
2839 template <int dim>
2841  const typename RefinementPossibilities<dim>::Possibilities refinement_case)
2842  : value(refinement_case)
2843 {
2844  // check that only those bits of
2845  // the given argument are set that
2846  // make sense for a given space
2847  // dimension
2848  Assert((refinement_case &
2850  refinement_case,
2851  ExcInvalidRefinementCase(refinement_case));
2852 }
2853 
2854 
2855 
2856 template <int dim>
2857 inline RefinementCase<dim>::RefinementCase(const std::uint8_t refinement_case)
2858  : value(refinement_case)
2859 {
2860  // check that only those bits of
2861  // the given argument are set that
2862  // make sense for a given space
2863  // dimension
2864  Assert((refinement_case &
2866  refinement_case,
2867  ExcInvalidRefinementCase(refinement_case));
2868 }
2869 
2870 
2871 
2872 template <int dim>
2873 inline DEAL_II_HOST_DEVICE RefinementCase<dim>::operator std::uint8_t() const
2874 {
2875  return value;
2876 }
2877 
2878 
2879 
2880 template <int dim>
2881 inline RefinementCase<dim>
2883 {
2884  return RefinementCase<dim>(static_cast<std::uint8_t>(value | r.value));
2885 }
2886 
2887 
2888 
2889 template <int dim>
2890 inline RefinementCase<dim>
2892 {
2893  return RefinementCase<dim>(static_cast<std::uint8_t>(value & r.value));
2894 }
2895 
2896 
2897 
2898 template <int dim>
2899 inline RefinementCase<dim>
2901 {
2902  return RefinementCase<dim>(static_cast<std::uint8_t>(
2904 }
2905 
2906 
2907 
2908 template <int dim>
2909 inline std::size_t
2911 {
2912  return sizeof(RefinementCase<dim>);
2913 }
2914 
2915 
2916 
2917 template <int dim>
2918 template <class Archive>
2919 inline void
2920 RefinementCase<dim>::serialize(Archive &ar, const unsigned int)
2921 {
2922  // serialization can't deal with bitfields, so copy from/to a full sized
2923  // std::uint8_t
2924  std::uint8_t uchar_value = value;
2925  ar &uchar_value;
2926  value = uchar_value;
2927 }
2928 
2929 
2930 
2931 template <>
2932 inline Point<1>
2933 GeometryInfo<1>::unit_cell_vertex(const unsigned int vertex)
2934 {
2935  AssertIndexRange(vertex, vertices_per_cell);
2936 
2937  return Point<1>(static_cast<double>(vertex));
2938 }
2939 
2940 
2941 
2942 template <>
2943 inline Point<2>
2944 GeometryInfo<2>::unit_cell_vertex(const unsigned int vertex)
2945 {
2946  AssertIndexRange(vertex, vertices_per_cell);
2947 
2948  return {static_cast<double>(vertex % 2), static_cast<double>(vertex / 2)};
2949 }
2950 
2951 
2952 
2953 template <>
2954 inline Point<3>
2955 GeometryInfo<3>::unit_cell_vertex(const unsigned int vertex)
2956 {
2957  AssertIndexRange(vertex, vertices_per_cell);
2958 
2959  return {static_cast<double>(vertex % 2),
2960  static_cast<double>(vertex / 2 % 2),
2961  static_cast<double>(vertex / 4)};
2962 }
2963 
2964 
2965 
2966 inline std::array<unsigned int, 0>
2968 {
2969  return {{}};
2970 }
2971 
2972 
2973 
2974 inline std::array<unsigned int, 1>
2976 {
2977  return {{0}};
2978 }
2979 
2980 
2981 
2982 template <int dim>
2985 {
2986  return {0U, faces_per_cell};
2987 }
2988 
2989 
2990 
2991 template <int dim>
2994 {
2995  return {0U, vertices_per_cell};
2996 }
2997 
2998 
2999 
3000 template <int dim>
3001 inline Point<dim>
3002 GeometryInfo<dim>::unit_cell_vertex(const unsigned int)
3003 {
3004  Assert(false, ExcNotImplemented());
3005 
3006  return {};
3007 }
3008 
3009 
3010 
3011 template <>
3012 inline unsigned int
3014 {
3015  Assert((p[0] >= 0) && (p[0] <= 1), ExcInvalidCoordinate(p[0]));
3016 
3017  return (p[0] <= 0.5 ? 0 : 1);
3018 }
3019 
3020 
3021 
3022 template <>
3023 inline unsigned int
3025 {
3026  Assert((p[0] >= 0) && (p[0] <= 1), ExcInvalidCoordinate(p[0]));
3027  Assert((p[1] >= 0) && (p[1] <= 1), ExcInvalidCoordinate(p[1]));
3028 
3029  return (p[0] <= 0.5 ? (p[1] <= 0.5 ? 0 : 2) : (p[1] <= 0.5 ? 1 : 3));
3030 }
3031 
3032 
3033 
3034 template <>
3035 inline unsigned int
3037 {
3038  Assert((p[0] >= 0) && (p[0] <= 1), ExcInvalidCoordinate(p[0]));
3039  Assert((p[1] >= 0) && (p[1] <= 1), ExcInvalidCoordinate(p[1]));
3040  Assert((p[2] >= 0) && (p[2] <= 1), ExcInvalidCoordinate(p[2]));
3041 
3042  return (p[0] <= 0.5 ?
3043  (p[1] <= 0.5 ? (p[2] <= 0.5 ? 0 : 4) : (p[2] <= 0.5 ? 2 : 6)) :
3044  (p[1] <= 0.5 ? (p[2] <= 0.5 ? 1 : 5) : (p[2] <= 0.5 ? 3 : 7)));
3045 }
3046 
3047 
3048 template <int dim>
3049 inline unsigned int
3051 {
3052  Assert(false, ExcNotImplemented());
3053 
3054  return 0;
3055 }
3056 
3057 
3058 
3059 template <>
3060 inline Point<1>
3062  const unsigned int child_index,
3063  const RefinementCase<1> refine_case)
3064 
3065 {
3066  AssertIndexRange(child_index, 2);
3068  (void)refine_case; // removes -Wunused-parameter warning in optimized mode
3069 
3070  return Point<1>(p * 2.0 - unit_cell_vertex(child_index));
3071 }
3072 
3073 
3074 
3075 template <>
3076 inline Point<2>
3078  const unsigned int child_index,
3079  const RefinementCase<2> refine_case)
3080 
3081 {
3082  AssertIndexRange(child_index, GeometryInfo<2>::n_children(refine_case));
3083 
3084  Point<2> point = p;
3085  switch (refine_case)
3086  {
3088  point[0] *= 2.0;
3089  if (child_index == 1)
3090  point[0] -= 1.0;
3091  break;
3093  point[1] *= 2.0;
3094  if (child_index == 1)
3095  point[1] -= 1.0;
3096  break;
3098  point *= 2.0;
3099  point -= unit_cell_vertex(child_index);
3100  break;
3101  default:
3102  Assert(false, ExcInternalError());
3103  }
3104 
3105  return point;
3106 }
3107 
3108 
3109 
3110 template <>
3111 inline Point<3>
3113  const unsigned int child_index,
3114  const RefinementCase<3> refine_case)
3115 
3116 {
3117  AssertIndexRange(child_index, GeometryInfo<3>::n_children(refine_case));
3118 
3119  Point<3> point = p;
3120  // there might be a cleverer way to do
3121  // this, but since this function is called
3122  // in very few places for initialization
3123  // purposes only, I don't care at the
3124  // moment
3125  switch (refine_case)
3126  {
3128  point[0] *= 2.0;
3129  if (child_index == 1)
3130  point[0] -= 1.0;
3131  break;
3133  point[1] *= 2.0;
3134  if (child_index == 1)
3135  point[1] -= 1.0;
3136  break;
3138  point[2] *= 2.0;
3139  if (child_index == 1)
3140  point[2] -= 1.0;
3141  break;
3143  point[0] *= 2.0;
3144  point[1] *= 2.0;
3145  if (child_index % 2 == 1)
3146  point[0] -= 1.0;
3147  if (child_index / 2 == 1)
3148  point[1] -= 1.0;
3149  break;
3151  // careful, this is slightly
3152  // different from xy and yz due to
3153  // different internal numbering of
3154  // children!
3155  point[0] *= 2.0;
3156  point[2] *= 2.0;
3157  if (child_index / 2 == 1)
3158  point[0] -= 1.0;
3159  if (child_index % 2 == 1)
3160  point[2] -= 1.0;
3161  break;
3163  point[1] *= 2.0;
3164  point[2] *= 2.0;
3165  if (child_index % 2 == 1)
3166  point[1] -= 1.0;
3167  if (child_index / 2 == 1)
3168  point[2] -= 1.0;
3169  break;
3171  point *= 2.0;
3172  point -= unit_cell_vertex(child_index);
3173  break;
3174  default:
3175  Assert(false, ExcInternalError());
3176  }
3177 
3178  return point;
3179 }
3180 
3181 
3182 
3183 template <int dim>
3184 inline Point<dim>
3186  const Point<dim> & /*p*/,
3187  const unsigned int /*child_index*/,
3188  const RefinementCase<dim> /*refine_case*/)
3189 
3190 {
3191  Assert(false, ExcNotImplemented());
3192  return {};
3193 }
3194 
3195 
3196 
3197 template <>
3198 inline Point<1>
3200  const unsigned int child_index,
3201  const RefinementCase<1> refine_case)
3202 
3203 {
3204  AssertIndexRange(child_index, 2);
3206  (void)refine_case; // removes -Wunused-parameter warning in optimized mode
3207 
3208  return (p + unit_cell_vertex(child_index)) * 0.5;
3209 }
3210 
3211 
3212 
3213 template <>
3214 inline Point<3>
3216  const unsigned int child_index,
3217  const RefinementCase<3> refine_case)
3218 
3219 {
3220  AssertIndexRange(child_index, GeometryInfo<3>::n_children(refine_case));
3221 
3222  Point<3> point = p;
3223  // there might be a cleverer way to do
3224  // this, but since this function is called
3225  // in very few places for initialization
3226  // purposes only, I don't care at the
3227  // moment
3228  switch (refine_case)
3229  {
3231  if (child_index == 1)
3232  point[0] += 1.0;
3233  point[0] *= 0.5;
3234  break;
3236  if (child_index == 1)
3237  point[1] += 1.0;
3238  point[1] *= 0.5;
3239  break;
3241  if (child_index == 1)
3242  point[2] += 1.0;
3243  point[2] *= 0.5;
3244  break;
3246  if (child_index % 2 == 1)
3247  point[0] += 1.0;
3248  if (child_index / 2 == 1)
3249  point[1] += 1.0;
3250  point[0] *= 0.5;
3251  point[1] *= 0.5;
3252  break;
3254  // careful, this is slightly
3255  // different from xy and yz due to
3256  // different internal numbering of
3257  // children!
3258  if (child_index / 2 == 1)
3259  point[0] += 1.0;
3260  if (child_index % 2 == 1)
3261  point[2] += 1.0;
3262  point[0] *= 0.5;
3263  point[2] *= 0.5;
3264  break;
3266  if (child_index % 2 == 1)
3267  point[1] += 1.0;
3268  if (child_index / 2 == 1)
3269  point[2] += 1.0;
3270  point[1] *= 0.5;
3271  point[2] *= 0.5;
3272  break;
3274  point += unit_cell_vertex(child_index);
3275  point *= 0.5;
3276  break;
3277  default:
3278  Assert(false, ExcInternalError());
3279  }
3280 
3281  return point;
3282 }
3283 
3284 
3285 
3286 template <>
3287 inline Point<2>
3289  const unsigned int child_index,
3290  const RefinementCase<2> refine_case)
3291 {
3292  AssertIndexRange(child_index, GeometryInfo<2>::n_children(refine_case));
3293 
3294  Point<2> point = p;
3295  switch (refine_case)
3296  {
3298  if (child_index == 1)
3299  point[0] += 1.0;
3300  point[0] *= 0.5;
3301  break;
3303  if (child_index == 1)
3304  point[1] += 1.0;
3305  point[1] *= 0.5;
3306  break;
3308  point += unit_cell_vertex(child_index);
3309  point *= 0.5;
3310  break;
3311  default:
3312  Assert(false, ExcInternalError());
3313  }
3314 
3315  return point;
3316 }
3317 
3318 
3319 
3320 template <int dim>
3321 inline Point<dim>
3323  const Point<dim> & /*p*/,
3324  const unsigned int /*child_index*/,
3325  const RefinementCase<dim> /*refine_case*/)
3326 {
3327  Assert(false, ExcNotImplemented());
3328  return {};
3329 }
3330 
3331 
3332 
3333 template <int dim>
3334 inline bool
3336 {
3337  Assert(false, ExcNotImplemented());
3338  return false;
3339 }
3340 
3341 template <>
3342 inline bool
3344 {
3345  return (p[0] >= 0.) && (p[0] <= 1.);
3346 }
3347 
3348 
3349 
3350 template <>
3351 inline bool
3353 {
3354  return (p[0] >= 0.) && (p[0] <= 1.) && (p[1] >= 0.) && (p[1] <= 1.);
3355 }
3356 
3357 
3358 
3359 template <>
3360 inline bool
3362 {
3363  return (p[0] >= 0.) && (p[0] <= 1.) && (p[1] >= 0.) && (p[1] <= 1.) &&
3364  (p[2] >= 0.) && (p[2] <= 1.);
3365 }
3366 
3367 
3368 
3369 template <int dim>
3370 inline bool
3372 {
3373  Assert(false, ExcNotImplemented());
3374  return false;
3375 }
3376 
3377 template <>
3378 inline bool
3379 GeometryInfo<1>::is_inside_unit_cell(const Point<1> &p, const double eps)
3380 {
3381  return (p[0] >= -eps) && (p[0] <= 1. + eps);
3382 }
3383 
3384 
3385 
3386 template <>
3387 inline bool
3388 GeometryInfo<2>::is_inside_unit_cell(const Point<2> &p, const double eps)
3389 {
3390  const double l = -eps, u = 1 + eps;
3391  return (p[0] >= l) && (p[0] <= u) && (p[1] >= l) && (p[1] <= u);
3392 }
3393 
3394 
3395 
3396 template <>
3397 inline bool
3398 GeometryInfo<3>::is_inside_unit_cell(const Point<3> &p, const double eps)
3399 {
3400  const double l = -eps, u = 1.0 + eps;
3401  return (p[0] >= l) && (p[0] <= u) && (p[1] >= l) && (p[1] <= u) &&
3402  (p[2] >= l) && (p[2] <= u);
3403 }
3404 
3405 
3406 
3407 template <>
3408 inline unsigned int
3409 GeometryInfo<1>::line_to_cell_vertices(const unsigned int line,
3410  const unsigned int vertex)
3411 {
3412  (void)line;
3413  AssertIndexRange(line, lines_per_cell);
3414  AssertIndexRange(vertex, 2);
3415 
3416  return vertex;
3417 }
3418 
3419 
3420 template <>
3421 inline unsigned int
3422 GeometryInfo<2>::line_to_cell_vertices(const unsigned int line,
3423  const unsigned int vertex)
3424 {
3425  constexpr unsigned int cell_vertices[4][2] = {{0, 2}, {1, 3}, {0, 1}, {2, 3}};
3426  return cell_vertices[line][vertex];
3427 }
3428 
3429 
3430 
3431 template <>
3432 inline unsigned int
3433 GeometryInfo<3>::line_to_cell_vertices(const unsigned int line,
3434  const unsigned int vertex)
3435 {
3436  AssertIndexRange(line, lines_per_cell);
3437  AssertIndexRange(vertex, 2);
3438 
3439  constexpr unsigned vertices[lines_per_cell][2] = {{0, 2}, // bottom face
3440  {1, 3},
3441  {0, 1},
3442  {2, 3},
3443  {4, 6}, // top face
3444  {5, 7},
3445  {4, 5},
3446  {6, 7},
3447  {0,
3448  4}, // connects of bottom
3449  {1, 5}, // top face
3450  {2, 6},
3451  {3, 7}};
3452  return vertices[line][vertex];
3453 }
3454 
3455 
3456 template <>
3457 inline unsigned int
3458 GeometryInfo<4>::line_to_cell_vertices(const unsigned int, const unsigned int)
3459 {
3460  Assert(false, ExcNotImplemented());
3462 }
3463 
3464 template <>
3465 inline unsigned int
3466 GeometryInfo<3>::standard_to_real_face_vertex(const unsigned int vertex,
3467  const bool face_orientation,
3468  const bool face_flip,
3469  const bool face_rotation)
3470 {
3472 
3473  // set up a table to make sure that
3474  // we handle non-standard faces correctly
3475  //
3476  // so set up a table that for each vertex (of
3477  // a quad in standard position) describes
3478  // which vertex to take
3479  //
3480  // first index: four vertices 0...3
3481  //
3482  // second index: face_orientation; 0:
3483  // opposite normal, 1: standard
3484  //
3485  // third index: face_flip; 0: standard, 1:
3486  // face rotated by 180 degrees
3487  //
3488  // forth index: face_rotation: 0: standard,
3489  // 1: face rotated by 90 degrees
3490 
3491  constexpr unsigned int vertex_translation[4][2][2][2] = {
3492  {{{0, 2}, // vertex 0, face_orientation=false, face_flip=false,
3493  // face_rotation=false and true
3494  {3, 1}}, // vertex 0, face_orientation=false, face_flip=true,
3495  // face_rotation=false and true
3496  {{0, 2}, // vertex 0, face_orientation=true, face_flip=false,
3497  // face_rotation=false and true
3498  {3, 1}}}, // vertex 0, face_orientation=true, face_flip=true,
3499  // face_rotation=false and true
3500 
3501  {{{2, 3}, // vertex 1 ...
3502  {1, 0}},
3503  {{1, 0}, {2, 3}}},
3504 
3505  {{{1, 0}, // vertex 2 ...
3506  {2, 3}},
3507  {{2, 3}, {1, 0}}},
3508 
3509  {{{3, 1}, // vertex 3 ...
3510  {0, 2}},
3511  {{3, 1}, {0, 2}}}};
3512 
3513  return vertex_translation[vertex][face_orientation][face_flip][face_rotation];
3514 }
3515 
3516 
3517 
3518 template <int dim>
3519 inline unsigned int
3520 GeometryInfo<dim>::standard_to_real_face_vertex(const unsigned int vertex,
3521  const bool,
3522  const bool,
3523  const bool)
3524 {
3525  Assert(dim > 1, ExcImpossibleInDim(dim));
3527  return vertex;
3528 }
3529 
3530 template <int dim>
3531 inline unsigned int
3533 {
3534  constexpr unsigned int n_children[RefinementCase<3>::cut_xyz + 1] = {
3535  0, 2, 2, 4, 2, 4, 4, 8};
3536 
3537  return n_children[ref_case];
3538 }
3539 
3540 
3541 
3542 template <int dim>
3543 inline unsigned int
3545 {
3546  Assert(false, ExcNotImplemented());
3547  return 0;
3548 }
3549 
3550 template <>
3551 inline unsigned int
3553 {
3554  Assert(false, ExcImpossibleInDim(1));
3555  return 0;
3556 }
3557 
3558 template <>
3559 inline unsigned int
3561 {
3562  return (subface_case == internal::SubfaceCase<2>::case_x) ? 2 : 0;
3563 }
3564 
3565 
3566 
3567 template <>
3568 inline unsigned int
3570 {
3571  const unsigned int nsubs[internal::SubfaceCase<3>::case_isotropic + 1] = {
3572  0, 2, 3, 3, 4, 2, 3, 3, 4, 4};
3573  return nsubs[subface_case];
3574 }
3575 
3576 
3577 
3578 template <int dim>
3579 inline double
3581  const unsigned int)
3582 {
3583  Assert(false, ExcNotImplemented());
3584  return 0.;
3585 }
3586 
3587 template <>
3588 inline double
3590  const unsigned int)
3591 {
3592  return 1;
3593 }
3594 
3595 
3596 template <>
3597 inline double
3599  const unsigned int)
3600 {
3601  double ratio = 1;
3602  switch (subface_case)
3603  {
3605  // Here, an
3606  // Assert(false,ExcInternalError())
3607  // would be the right
3608  // choice, but
3609  // unfortunately the
3610  // current function is
3611  // also called for faces
3612  // without children (see
3613  // tests/fe/mapping.cc).
3614  // Assert(false, ExcMessage("Face has no subfaces."));
3615  // Furthermore, assign
3616  // following value as
3617  // otherwise the
3618  // bits/volume_x tests
3619  // break
3621  break;
3623  ratio = 0.5;
3624  break;
3625  default:
3626  // there should be no
3627  // cases left
3628  Assert(false, ExcInternalError());
3629  break;
3630  }
3631 
3632  return ratio;
3633 }
3634 
3635 
3636 template <>
3637 inline double
3639  const unsigned int subface_no)
3640 {
3641  double ratio = 1;
3642  switch (subface_case)
3643  {
3645  // Here, an
3646  // Assert(false,ExcInternalError())
3647  // would be the right
3648  // choice, but
3649  // unfortunately the
3650  // current function is
3651  // also called for faces
3652  // without children (see
3653  // tests/bits/mesh_3d_16.cc). Add
3654  // following switch to
3655  // avoid diffs in
3656  // tests/bits/mesh_3d_16
3658  break;
3661  ratio = 0.5;
3662  break;
3666  ratio = 0.25;
3667  break;
3670  if (subface_no < 2)
3671  ratio = 0.25;
3672  else
3673  ratio = 0.5;
3674  break;
3677  if (subface_no == 0)
3678  ratio = 0.5;
3679  else
3680  ratio = 0.25;
3681  break;
3682  default:
3683  // there should be no
3684  // cases left
3685  Assert(false, ExcInternalError());
3686  break;
3687  }
3688 
3689  return ratio;
3690 }
3691 
3692 
3693 
3694 template <int dim>
3696  const RefinementCase<dim> &,
3697  const unsigned int,
3698  const bool,
3699  const bool,
3700  const bool)
3701 {
3702  Assert(false, ExcNotImplemented());
3703  return RefinementCase<dim - 1>::no_refinement;
3704 }
3705 
3706 template <>
3708  const RefinementCase<1> &,
3709  const unsigned int,
3710  const bool,
3711  const bool,
3712  const bool)
3713 {
3714  Assert(false, ExcImpossibleInDim(1));
3715 
3717 }
3718 
3719 
3720 template <>
3721 inline RefinementCase<1>
3723  const RefinementCase<2> &cell_refinement_case,
3724  const unsigned int face_no,
3725  const bool,
3726  const bool,
3727  const bool)
3728 {
3729  const unsigned int dim = 2;
3730  AssertIndexRange(cell_refinement_case,
3733 
3734  // simple special case
3735  if (cell_refinement_case == RefinementCase<dim>::cut_xy)
3736  return RefinementCase<1>::cut_x;
3737 
3738  const RefinementCase<dim - 1>
3741  {RefinementCase<dim - 1>::no_refinement, // no_refinement
3742  RefinementCase<dim - 1>::no_refinement},
3743 
3744  {RefinementCase<dim - 1>::no_refinement, RefinementCase<dim - 1>::cut_x},
3745 
3746  {RefinementCase<dim - 1>::cut_x, RefinementCase<dim - 1>::no_refinement},
3747 
3748  {RefinementCase<dim - 1>::cut_x, // cut_xy
3749  RefinementCase<dim - 1>::cut_x}};
3750 
3751  return ref_cases[cell_refinement_case][face_no / 2];
3752 }
3753 
3754 
3755 template <>
3756 inline RefinementCase<2>
3758  const RefinementCase<3> &cell_refinement_case,
3759  const unsigned int face_no,
3760  const bool face_orientation,
3761  const bool /*face_flip*/,
3762  const bool face_rotation)
3763 {
3764  const unsigned int dim = 3;
3765  AssertIndexRange(cell_refinement_case,
3768 
3769  // simple special case
3770  if (cell_refinement_case == RefinementCase<dim>::cut_xyz)
3771  return RefinementCase<dim - 1>::cut_xy;
3772 
3773  const RefinementCase<dim - 1>
3776  {RefinementCase<dim - 1>::no_refinement, // no_refinement
3777  RefinementCase<dim - 1>::no_refinement,
3778  RefinementCase<dim - 1>::no_refinement},
3779 
3780  {RefinementCase<dim - 1>::no_refinement, // cut_x
3781  RefinementCase<dim - 1>::cut_y,
3782  RefinementCase<dim - 1>::cut_x},
3783 
3784  {RefinementCase<dim - 1>::cut_x, // cut_y
3785  RefinementCase<dim - 1>::no_refinement,
3786  RefinementCase<dim - 1>::cut_y},
3787 
3788  {RefinementCase<dim - 1>::cut_x, // cut_xy
3789  RefinementCase<dim - 1>::cut_y,
3790  RefinementCase<dim - 1>::cut_xy},
3791 
3792  {RefinementCase<dim - 1>::cut_y, // cut_z
3793  RefinementCase<dim - 1>::cut_x,
3794  RefinementCase<dim - 1>::no_refinement},
3795 
3796  {RefinementCase<dim - 1>::cut_y, // cut_xz
3797  RefinementCase<dim - 1>::cut_xy,
3798  RefinementCase<dim - 1>::cut_x},
3799 
3800  {RefinementCase<dim - 1>::cut_xy, // cut_yz
3801  RefinementCase<dim - 1>::cut_x,
3802  RefinementCase<dim - 1>::cut_y},
3803 
3804  {RefinementCase<dim - 1>::cut_xy, // cut_xyz
3805  RefinementCase<dim - 1>::cut_xy,
3806  RefinementCase<dim - 1>::cut_xy},
3807  };
3808 
3809  const RefinementCase<dim - 1> ref_case =
3810  ref_cases[cell_refinement_case][face_no / 2];
3811 
3812  const RefinementCase<dim - 1> flip[4] = {
3813  RefinementCase<dim - 1>::no_refinement,
3814  RefinementCase<dim - 1>::cut_y,
3815  RefinementCase<dim - 1>::cut_x,
3816  RefinementCase<dim - 1>::cut_xy};
3817 
3818  // correct the ref_case for face_orientation
3819  // and face_rotation. for face_orientation,
3820  // 'true' is the default value whereas for
3821  // face_rotation, 'false' is standard. If
3822  // <tt>face_rotation==face_orientation</tt>,
3823  // then one of them is non-standard and we
3824  // have to swap cut_x and cut_y, otherwise no
3825  // change is necessary. face_flip has no
3826  // influence. however, in order to keep the
3827  // interface consistent with other functions,
3828  // we still include it as an argument to this
3829  // function
3830  return (face_orientation == face_rotation) ? flip[ref_case] : ref_case;
3831 }
3832 
3833 
3834 
3835 template <int dim>
3836 inline RefinementCase<1>
3838  const unsigned int)
3839 {
3840  Assert(false, ExcNotImplemented());
3842 }
3843 
3844 template <>
3845 inline RefinementCase<1>
3847  const RefinementCase<1> &cell_refinement_case,
3848  const unsigned int line_no)
3849 {
3850  (void)line_no;
3851  const unsigned int dim = 1;
3852  (void)dim;
3853  AssertIndexRange(cell_refinement_case,
3856 
3857  return cell_refinement_case;
3858 }
3859 
3860 
3861 template <>
3862 inline RefinementCase<1>
3864  const RefinementCase<2> &cell_refinement_case,
3865  const unsigned int line_no)
3866 {
3867  // Assertions are in face_refinement_case()
3868  return face_refinement_case(cell_refinement_case, line_no);
3869 }
3870 
3871 
3872 template <>
3873 inline RefinementCase<1>
3875  const RefinementCase<3> &cell_refinement_case,
3876  const unsigned int line_no)
3877 {
3878  const unsigned int dim = 3;
3879  AssertIndexRange(cell_refinement_case,
3882 
3883  // simple special case
3884  if (cell_refinement_case == RefinementCase<dim>::cut_xyz)
3885  return RefinementCase<1>::cut_x;
3886 
3887  // array indicating, which simple refine
3888  // case cuts a line in direction x, y or
3889  // z. For example, cut_y and everything
3890  // containing cut_y (cut_xy, cut_yz,
3891  // cut_xyz) cuts lines, which are in y
3892  // direction.
3893  const RefinementCase<dim> cut_one[dim] = {RefinementCase<dim>::cut_x,
3896 
3897  // order the direction of lines
3898  // 0->x, 1->y, 2->z
3899  const unsigned int direction[lines_per_cell] = {
3900  1, 1, 0, 0, 1, 1, 0, 0, 2, 2, 2, 2};
3901 
3902  return ((cell_refinement_case & cut_one[direction[line_no]]) ?
3905 }
3906 
3907 
3908 
3909 template <int dim>
3910 inline RefinementCase<dim>
3912  const RefinementCase<dim - 1> &,
3913  const unsigned int,
3914  const bool,
3915  const bool,
3916  const bool)
3917 {
3918  Assert(false, ExcNotImplemented());
3919 
3921 }
3922 
3923 template <>
3924 inline RefinementCase<1>
3926  const RefinementCase<0> &,
3927  const unsigned int,
3928  const bool,
3929  const bool,
3930  const bool)
3931 {
3932  const unsigned int dim = 1;
3933  Assert(false, ExcImpossibleInDim(dim));
3934 
3936 }
3937 
3938 
3939 template <>
3940 inline RefinementCase<2>
3942  const RefinementCase<1> &face_refinement_case,
3943  const unsigned int face_no,
3944  const bool,
3945  const bool,
3946  const bool)
3947 {
3948  const unsigned int dim = 2;
3949  AssertIndexRange(face_refinement_case,
3952 
3953  if (face_refinement_case == RefinementCase<dim>::cut_x)
3954  return (face_no / 2) != 0u ? RefinementCase<dim>::cut_x :
3956  else
3958 }
3959 
3960 
3961 template <>
3962 inline RefinementCase<3>
3964  const RefinementCase<2> &face_refinement_case,
3965  const unsigned int face_no,
3966  const bool face_orientation,
3967  const bool /*face_flip*/,
3968  const bool face_rotation)
3969 {
3970  const unsigned int dim = 3;
3971  AssertIndexRange(face_refinement_case,
3974 
3979 
3980  // correct the face_refinement_case for
3981  // face_orientation and face_rotation. for
3982  // face_orientation, 'true' is the default
3983  // value whereas for face_rotation, 'false'
3984  // is standard. If
3985  // <tt>face_rotation==face_orientation</tt>,
3986  // then one of them is non-standard and we
3987  // have to swap cut_x and cut_y, otherwise no
3988  // change is necessary. face_flip has no
3989  // influence. however, in order to keep the
3990  // interface consistent with other functions,
3991  // we still include it as an argument to this
3992  // function
3993  const RefinementCase<dim - 1> std_face_ref =
3994  (face_orientation == face_rotation) ? flip[face_refinement_case] :
3995  face_refinement_case;
3996 
3997  const RefinementCase<dim> face_to_cell[3][4] = {
3998  {RefinementCase<dim>::no_refinement, // faces 0 and 1
3999  RefinementCase<dim>::cut_y, // cut_x in face 0 means cut_y for the cell
4002 
4003  {RefinementCase<dim>::no_refinement, // faces 2 and 3 (note that x and y are
4004  // "exchanged on faces 2 and 3")
4008 
4009  {RefinementCase<dim>::no_refinement, // faces 4 and 5
4013 
4014  return face_to_cell[face_no / 2][std_face_ref];
4015 }
4016 
4017 
4018 
4019 template <int dim>
4020 inline RefinementCase<dim>
4022  const unsigned int)
4023 {
4024  Assert(false, ExcNotImplemented());
4025 
4027 }
4028 
4029 template <>
4030 inline RefinementCase<1>
4032  const unsigned int line_no)
4033 {
4034  (void)line_no;
4035  AssertIndexRange(line_no, 1);
4036 
4037  return RefinementCase<1>::cut_x;
4038 }
4039 
4040 
4041 template <>
4042 inline RefinementCase<2>
4044  const unsigned int line_no)
4045 {
4046  const unsigned int dim = 2;
4047  (void)dim;
4049 
4050  return (line_no / 2) != 0u ? RefinementCase<2>::cut_x :
4052 }
4053 
4054 
4055 template <>
4056 inline RefinementCase<3>
4058  const unsigned int line_no)
4059 {
4060  const unsigned int dim = 3;
4062 
4063  const RefinementCase<dim> ref_cases[6] = {
4064  RefinementCase<dim>::cut_y, // lines 0 and 1
4065  RefinementCase<dim>::cut_x, // lines 2 and 3
4066  RefinementCase<dim>::cut_y, // lines 4 and 5
4067  RefinementCase<dim>::cut_x, // lines 6 and 7
4068  RefinementCase<dim>::cut_z, // lines 8 and 9
4069  RefinementCase<dim>::cut_z}; // lines 10 and 11
4070 
4071  return ref_cases[line_no / 2];
4072 }
4073 
4074 
4075 
4076 template <>
4077 inline unsigned int
4078 GeometryInfo<3>::real_to_standard_face_vertex(const unsigned int vertex,
4079  const bool face_orientation,
4080  const bool face_flip,
4081  const bool face_rotation)
4082 {
4084 
4085  // set up a table to make sure that
4086  // we handle non-standard faces correctly
4087  //
4088  // so set up a table that for each vertex (of
4089  // a quad in standard position) describes
4090  // which vertex to take
4091  //
4092  // first index: four vertices 0...3
4093  //
4094  // second index: face_orientation; 0:
4095  // opposite normal, 1: standard
4096  //
4097  // third index: face_flip; 0: standard, 1:
4098  // face rotated by 180 degrees
4099  //
4100  // forth index: face_rotation: 0: standard,
4101  // 1: face rotated by 90 degrees
4102 
4103  const unsigned int vertex_translation[4][2][2][2] = {
4104  {{{0, 2}, // vertex 0, face_orientation=false, face_flip=false,
4105  // face_rotation=false and true
4106  {3, 1}}, // vertex 0, face_orientation=false, face_flip=true,
4107  // face_rotation=false and true
4108  {{0, 1}, // vertex 0, face_orientation=true, face_flip=false,
4109  // face_rotation=false and true
4110  {3, 2}}}, // vertex 0, face_orientation=true, face_flip=true,
4111  // face_rotation=false and true
4112 
4113  {{{2, 3}, // vertex 1 ...
4114  {1, 0}},
4115  {{1, 3}, {2, 0}}},
4116 
4117  {{{1, 0}, // vertex 2 ...
4118  {2, 3}},
4119  {{2, 0}, {1, 3}}},
4120 
4121  {{{3, 1}, // vertex 3 ...
4122  {0, 2}},
4123  {{3, 2}, {0, 1}}}};
4124 
4125  return vertex_translation[vertex][face_orientation][face_flip][face_rotation];
4126 }
4127 
4128 
4129 
4130 template <int dim>
4131 inline unsigned int
4132 GeometryInfo<dim>::real_to_standard_face_vertex(const unsigned int vertex,
4133  const bool,
4134  const bool,
4135  const bool)
4136 {
4137  Assert(dim > 1, ExcImpossibleInDim(dim));
4139  return vertex;
4140 }
4141 
4142 
4143 
4144 template <>
4145 inline unsigned int
4146 GeometryInfo<3>::standard_to_real_face_line(const unsigned int line,
4147  const bool face_orientation,
4148  const bool face_flip,
4149  const bool face_rotation)
4150 {
4152 
4153 
4154  // make sure we handle
4155  // non-standard faces correctly
4156  //
4157  // so set up a table that for each line (of a
4158  // quad) describes which line to take
4159  //
4160  // first index: four lines 0...3
4161  //
4162  // second index: face_orientation; 0:
4163  // opposite normal, 1: standard
4164  //
4165  // third index: face_flip; 0: standard, 1:
4166  // face rotated by 180 degrees
4167  //
4168  // forth index: face_rotation: 0: standard,
4169  // 1: face rotated by 90 degrees
4170 
4171  const unsigned int line_translation[4][2][2][2] = {
4172  {{{2, 0}, // line 0, face_orientation=false, face_flip=false,
4173  // face_rotation=false and true
4174  {3, 1}}, // line 0, face_orientation=false, face_flip=true,
4175  // face_rotation=false and true
4176  {{0, 3}, // line 0, face_orientation=true, face_flip=false,
4177  // face_rotation=false and true
4178  {1, 2}}}, // line 0, face_orientation=true, face_flip=true,
4179  // face_rotation=false and true
4180 
4181  {{{3, 1}, // line 1 ...
4182  {2, 0}},
4183  {{1, 2}, {0, 3}}},
4184 
4185  {{{0, 3}, // line 2 ...
4186  {1, 2}},
4187  {{2, 0}, {3, 1}}},
4188 
4189  {{{1, 2}, // line 3 ...
4190  {0, 3}},
4191  {{3, 1}, {2, 0}}}};
4192 
4193  return line_translation[line][face_orientation][face_flip][face_rotation];
4194 }
4195 
4196 
4197 
4198 template <int dim>
4199 inline unsigned int
4200 GeometryInfo<dim>::standard_to_real_face_line(const unsigned int line,
4201  const bool,
4202  const bool,
4203  const bool)
4204 {
4205  Assert(false, ExcNotImplemented());
4206  return line;
4207 }
4208 
4209 
4210 
4211 template <>
4212 inline unsigned int
4213 GeometryInfo<2>::standard_to_real_line_vertex(const unsigned int vertex,
4214  const bool line_orientation)
4215 {
4216  return line_orientation ? vertex : (1 - vertex);
4217 }
4218 
4219 
4220 
4221 template <int dim>
4222 inline unsigned int
4223 GeometryInfo<dim>::standard_to_real_line_vertex(const unsigned int vertex,
4224  const bool)
4225 {
4226  Assert(false, ExcNotImplemented());
4227  return vertex;
4228 }
4229 
4230 
4231 
4232 template <>
4233 inline std::array<unsigned int, 2>
4235  const unsigned int vertex)
4236 {
4237  return {{vertex % 2, vertex / 2}};
4238 }
4239 
4240 
4241 
4242 template <int dim>
4243 inline std::array<unsigned int, 2>
4245  const unsigned int vertex)
4246 {
4247  Assert(false, ExcNotImplemented());
4248  (void)vertex;
4249  return {{0, 0}};
4250 }
4251 
4252 
4253 
4254 template <>
4255 inline std::array<unsigned int, 2>
4257 {
4258  // set up a table that for each
4259  // line describes a) from which
4260  // quad to take it, b) which line
4261  // therein it is if the face is
4262  // oriented correctly
4263  static const unsigned int lookup_table[GeometryInfo<3>::lines_per_cell][2] = {
4264  {4, 0}, // take first four lines from bottom face
4265  {4, 1},
4266  {4, 2},
4267  {4, 3},
4268 
4269  {5, 0}, // second four lines from top face
4270  {5, 1},
4271  {5, 2},
4272  {5, 3},
4273 
4274  {0, 0}, // the rest randomly
4275  {1, 0},
4276  {0, 1},
4277  {1, 1}};
4278 
4279  return {{lookup_table[i][0], lookup_table[i][1]}};
4280 }
4281 
4282 
4283 
4284 template <int dim>
4285 inline std::array<unsigned int, 2>
4287 {
4288  Assert(false, ExcNotImplemented());
4289  (void)line;
4290  return {{0, 0}};
4291 }
4292 
4293 
4294 
4295 template <>
4296 inline std::array<unsigned int, 2>
4298  const unsigned int vertex)
4299 {
4300  // get the corner indices by asking either the bottom or the top face for its
4301  // vertices. handle non-standard faces by calling the vertex reordering
4302  // function from GeometryInfo
4303 
4304  // bottom face (4) for first four vertices, top face (5) for the rest
4305  return {{4 + vertex / 4, vertex % 4}};
4306 }
4307 
4308 
4309 
4310 template <int dim>
4311 inline std::array<unsigned int, 2>
4313  const unsigned int vertex)
4314 {
4315  Assert(false, ExcNotImplemented());
4316  (void)vertex;
4317  return {{0, 0}};
4318 }
4319 
4320 
4321 
4322 template <>
4323 inline unsigned int
4324 GeometryInfo<3>::real_to_standard_face_line(const unsigned int line,
4325  const bool face_orientation,
4326  const bool face_flip,
4327  const bool face_rotation)
4328 {
4330 
4331 
4332  // make sure we handle
4333  // non-standard faces correctly
4334  //
4335  // so set up a table that for each line (of a
4336  // quad) describes which line to take
4337  //
4338  // first index: four lines 0...3
4339  //
4340  // second index: face_orientation; 0:
4341  // opposite normal, 1: standard
4342  //
4343  // third index: face_flip; 0: standard, 1:
4344  // face rotated by 180 degrees
4345  //
4346  // forth index: face_rotation: 0: standard,
4347  // 1: face rotated by 90 degrees
4348 
4349  const unsigned int line_translation[4][2][2][2] = {
4350  {{{2, 0}, // line 0, face_orientation=false, face_flip=false,
4351  // face_rotation=false and true
4352  {3, 1}}, // line 0, face_orientation=false, face_flip=true,
4353  // face_rotation=false and true
4354  {{0, 2}, // line 0, face_orientation=true, face_flip=false,
4355  // face_rotation=false and true
4356  {1, 3}}}, // line 0, face_orientation=true, face_flip=true,
4357  // face_rotation=false and true
4358 
4359  {{{3, 1}, // line 1 ...
4360  {2, 0}},
4361  {{1, 3}, {0, 2}}},
4362 
4363  {{{0, 3}, // line 2 ...
4364  {1, 2}},
4365  {{2, 1}, {3, 0}}},
4366 
4367  {{{1, 2}, // line 3 ...
4368  {0, 3}},
4369  {{3, 0}, {2, 1}}}};
4370 
4371  return line_translation[line][face_orientation][face_flip][face_rotation];
4372 }
4373 
4374 
4375 
4376 template <int dim>
4377 inline unsigned int
4378 GeometryInfo<dim>::real_to_standard_face_line(const unsigned int line,
4379  const bool,
4380  const bool,
4381  const bool)
4382 {
4383  Assert(false, ExcNotImplemented());
4384  return line;
4385 }
4386 
4387 
4388 
4389 template <>
4390 inline unsigned int
4392  const unsigned int face,
4393  const unsigned int subface,
4394  const bool,
4395  const bool,
4396  const bool,
4397  const RefinementCase<0> &)
4398 {
4399  (void)subface;
4400  AssertIndexRange(face, faces_per_cell);
4401  AssertIndexRange(subface, max_children_per_face);
4402 
4403  return face;
4404 }
4405 
4406 
4407 
4408 template <>
4409 inline unsigned int
4411  const unsigned int face,
4412  const unsigned int subface,
4413  const bool /*face_orientation*/,
4414  const bool face_flip,
4415  const bool /*face_rotation*/,
4416  const RefinementCase<1> &)
4417 {
4418  AssertIndexRange(face, faces_per_cell);
4419  AssertIndexRange(subface, max_children_per_face);
4420 
4421  // always return the child adjacent to the specified
4422  // subface. if the face of a cell is not refined, don't
4423  // throw an assertion but deliver the child adjacent to
4424  // the face nevertheless, i.e. deliver the child of
4425  // this cell adjacent to the subface of a possibly
4426  // refined neighbor. this simplifies setting neighbor
4427  // information in execute_refinement.
4428  constexpr unsigned int
4429  subcells[/* possible face orientation */ 2]
4430  [/* number of different ways to refine a cell */ 4]
4431  [/* faces_per_cell */ 4][/* max_children_per_face */ 2] = {
4432  {
4433  // Normal orientation (face_flip = false)
4434  {{0, 0}, {1, 1}, {0, 1}, {0, 1}}, // cut_x
4435  {{0, 1}, {0, 1}, {0, 0}, {1, 1}}, // cut_y
4436  {{0, 2}, {1, 3}, {0, 1}, {2, 3}} // cut_xy, i.e., isotropic
4437  },
4438  {
4439  // Flipped orientation (face_flip = true)
4440  {{0, 0}, {1, 1}, {1, 0}, {1, 0}}, // cut_x
4441  {{1, 0}, {1, 0}, {0, 0}, {1, 1}}, // cut_y
4442  {{2, 0}, {3, 1}, {1, 0}, {3, 2}} // cut_xy, i.e., isotropic
4443  }};
4444 
4445  return subcells[face_flip][ref_case - 1][face][subface];
4446 }
4447 
4448 
4449 
4450 template <>
4451 inline unsigned int
4453  const unsigned int face,
4454  const unsigned int subface,
4455  const bool face_orientation,
4456  const bool face_flip,
4457  const bool face_rotation,
4458  const RefinementCase<2> &face_ref_case)
4459 {
4460  const unsigned int dim = 3;
4461 
4463  ExcMessage("Cell has no children."));
4464  AssertIndexRange(face, faces_per_cell);
4465  if (!(subface == 0 &&
4466  face_ref_case == RefinementCase<dim - 1>::no_refinement))
4467  {
4468  AssertIndexRange(subface,
4469  GeometryInfo<dim - 1>::n_children(face_ref_case));
4470  }
4471 
4472  // invalid number used for invalid cases,
4473  // e.g. when the children are more refined at
4474  // a given face than the face itself
4475  const unsigned int e = numbers::invalid_unsigned_int;
4476 
4477  // the whole process of finding a child cell
4478  // at a given subface considering the
4479  // possibly anisotropic refinement cases of
4480  // the cell and the face as well as
4481  // orientation, flip and rotation of the face
4482  // is quite complicated. thus, we break it
4483  // down into several steps.
4484 
4485  // first step: convert the given face refine
4486  // case to a face refine case concerning the
4487  // face in standard orientation (, flip and
4488  // rotation). This only affects cut_x and
4489  // cut_y
4490  const RefinementCase<dim - 1> flip[4] = {
4491  RefinementCase<dim - 1>::no_refinement,
4492  RefinementCase<dim - 1>::cut_y,
4493  RefinementCase<dim - 1>::cut_x,
4494  RefinementCase<dim - 1>::cut_xy};
4495  // for face_orientation, 'true' is the
4496  // default value whereas for face_rotation,
4497  // 'false' is standard. If
4498  // <tt>face_rotation==face_orientation</tt>,
4499  // then one of them is non-standard and we
4500  // have to swap cut_x and cut_y, otherwise no
4501  // change is necessary.
4502  const RefinementCase<dim - 1> std_face_ref =
4503  (face_orientation == face_rotation) ? flip[face_ref_case] : face_ref_case;
4504 
4505  // second step: convert the given subface
4506  // index to the one for a standard face
4507  // respecting face_orientation, face_flip and
4508  // face_rotation
4509 
4510  // first index: face_ref_case
4511  // second index: face_orientation
4512  // third index: face_flip
4513  // forth index: face_rotation
4514  // fifth index: subface index
4515  const unsigned int subface_exchange[4][2][2][2][4] = {
4516  // no_refinement (subface 0 stays 0,
4517  // all others are invalid)
4518  {{{{0, e, e, e}, {0, e, e, e}}, {{0, e, e, e}, {0, e, e, e}}},
4519  {{{0, e, e, e}, {0, e, e, e}}, {{0, e, e, e}, {0, e, e, e}}}},
4520  // cut_x (here, if the face is only
4521  // rotated OR only falsely oriented,
4522  // then subface 0 of the non-standard
4523  // face does NOT correspond to one of
4524  // the subfaces of a standard
4525  // face. Thus we indicate the subface
4526  // which is located at the lower left
4527  // corner (the origin of the face's
4528  // local coordinate system) with
4529  // '0'. The rest of this issue is
4530  // taken care of using the above
4531  // conversion to a 'standard face
4532  // refine case')
4533  {{{{0, 1, e, e}, {0, 1, e, e}}, {{1, 0, e, e}, {1, 0, e, e}}},
4534  {{{0, 1, e, e}, {0, 1, e, e}}, {{1, 0, e, e}, {1, 0, e, e}}}},
4535  // cut_y (the same applies as for
4536  // cut_x)
4537  {{{{0, 1, e, e}, {1, 0, e, e}}, {{1, 0, e, e}, {0, 1, e, e}}},
4538  {{{0, 1, e, e}, {1, 0, e, e}}, {{1, 0, e, e}, {0, 1, e, e}}}},
4539  // cut_xyz: this information is
4540  // identical to the information
4541  // returned by
4542  // GeometryInfo<3>::real_to_standard_face_vertex()
4543  {{{{0, 2, 1, 3}, // face_orientation=false, face_flip=false,
4544  // face_rotation=false, subfaces 0,1,2,3
4545  {2, 3, 0, 1}}, // face_orientation=false, face_flip=false,
4546  // face_rotation=true, subfaces 0,1,2,3
4547  {{3, 1, 2, 0}, // face_orientation=false, face_flip=true,
4548  // face_rotation=false, subfaces 0,1,2,3
4549  {1, 0, 3, 2}}}, // face_orientation=false, face_flip=true,
4550  // face_rotation=true, subfaces 0,1,2,3
4551  {{{0, 1, 2, 3}, // face_orientation=true, face_flip=false,
4552  // face_rotation=false, subfaces 0,1,2,3
4553  {1, 3, 0, 2}}, // face_orientation=true, face_flip=false,
4554  // face_rotation=true, subfaces 0,1,2,3
4555  {{3, 2, 1, 0}, // face_orientation=true, face_flip=true,
4556  // face_rotation=false, subfaces 0,1,2,3
4557  {2, 0, 3, 1}}}}}; // face_orientation=true, face_flip=true,
4558  // face_rotation=true, subfaces 0,1,2,3
4559 
4560  const unsigned int std_subface =
4561  subface_exchange[face_ref_case][face_orientation][face_flip][face_rotation]
4562  [subface];
4563  Assert(std_subface != e, ExcInternalError());
4564 
4565  // third step: these are the children, which
4566  // can be found at the given subfaces of an
4567  // isotropically refined (standard) face
4568  //
4569  // first index: (refinement_case-1)
4570  // second index: face_index
4571  // third index: subface_index (isotropic refinement)
4572  const unsigned int iso_children[RefinementCase<dim>::cut_xyz][faces_per_cell]
4573  [max_children_per_face] = {
4574  // cut_x
4575  {{0, 0, 0, 0}, // face 0, subfaces 0,1,2,3
4576  {1, 1, 1, 1}, // face 1, subfaces 0,1,2,3
4577  {0, 0, 1, 1}, // face 2, subfaces 0,1,2,3
4578  {0, 0, 1, 1}, // face 3, subfaces 0,1,2,3
4579  {0, 1, 0, 1}, // face 4, subfaces 0,1,2,3
4580  {0, 1, 0, 1}}, // face 5, subfaces 0,1,2,3
4581  // cut_y
4582  {{0, 1, 0, 1},
4583  {0, 1, 0, 1},
4584  {0, 0, 0, 0},
4585  {1, 1, 1, 1},
4586  {0, 0, 1, 1},
4587  {0, 0, 1, 1}},
4588  // cut_xy
4589  {{0, 2, 0, 2},
4590  {1, 3, 1, 3},
4591  {0, 0, 1, 1},
4592  {2, 2, 3, 3},
4593  {0, 1, 2, 3},
4594  {0, 1, 2, 3}},
4595  // cut_z
4596  {{0, 0, 1, 1},
4597  {0, 0, 1, 1},
4598  {0, 1, 0, 1},
4599  {0, 1, 0, 1},
4600  {0, 0, 0, 0},
4601  {1, 1, 1, 1}},
4602  // cut_xz
4603  {{0, 0, 1, 1},
4604  {2, 2, 3, 3},
4605  {0, 1, 2, 3},
4606  {0, 1, 2, 3},
4607  {0, 2, 0, 2},
4608  {1, 3, 1, 3}},
4609  // cut_yz
4610  {{0, 1, 2, 3},
4611  {0, 1, 2, 3},
4612  {0, 2, 0, 2},
4613  {1, 3, 1, 3},
4614  {0, 0, 1, 1},
4615  {2, 2, 3, 3}},
4616  // cut_xyz
4617  {{0, 2, 4, 6},
4618  {1, 3, 5, 7},
4619  {0, 4, 1, 5},
4620  {2, 6, 3, 7},
4621  {0, 1, 2, 3},
4622  {4, 5, 6, 7}}};
4623 
4624  // forth step: check, whether the given face
4625  // refine case is valid for the given cell
4626  // refine case. this is the case, if the
4627  // given face refine case is at least as
4628  // refined as the face is for the given cell
4629  // refine case
4630 
4631  // note, that we are considering standard
4632  // face refinement cases here and thus must
4633  // not pass the given orientation, flip and
4634  // rotation flags
4635  if ((std_face_ref & face_refinement_case(ref_case, face)) ==
4636  face_refinement_case(ref_case, face))
4637  {
4638  // all is fine. for anisotropic face
4639  // refine cases, select one of the
4640  // isotropic subfaces which neighbors the
4641  // same child
4642 
4643  // first index: (standard) face refine case
4644  // second index: subface index
4645  const unsigned int equivalent_iso_subface[4][4] = {
4646  {0, e, e, e}, // no_refinement
4647  {0, 3, e, e}, // cut_x
4648  {0, 3, e, e}, // cut_y
4649  {0, 1, 2, 3}}; // cut_xy
4650 
4651  const unsigned int equ_std_subface =
4652  equivalent_iso_subface[std_face_ref][std_subface];
4653  Assert(equ_std_subface != e, ExcInternalError());
4654 
4655  return iso_children[ref_case - 1][face][equ_std_subface];
4656  }
4657  else
4658  {
4659  // the face_ref_case was too coarse,
4660  // throw an error
4661  Assert(false,
4662  ExcMessage("The face RefineCase is too coarse "
4663  "for the given cell RefineCase."));
4664  }
4665  // we only get here in case of an error
4666  return e;
4667 }
4668 
4669 
4670 
4671 template <>
4672 inline unsigned int
4674  const unsigned int,
4675  const unsigned int,
4676  const bool,
4677  const bool,
4678  const bool,
4679  const RefinementCase<3> &)
4680 {
4681  Assert(false, ExcNotImplemented());
4683 }
4684 
4685 
4686 
4687 template <>
4688 inline unsigned int
4689 GeometryInfo<2>::face_to_cell_lines(const unsigned int face,
4690  const unsigned int line,
4691  const bool,
4692  const bool,
4693  const bool)
4694 {
4695  (void)line;
4696  AssertIndexRange(face, faces_per_cell);
4697  AssertIndexRange(line, lines_per_face);
4698 
4699  // The face is a line itself.
4700  return face;
4701 }
4702 
4703 
4704 
4705 template <>
4706 inline unsigned int
4707 GeometryInfo<3>::face_to_cell_lines(const unsigned int face,
4708  const unsigned int line,
4709  const bool face_orientation,
4710  const bool face_flip,
4711  const bool face_rotation)
4712 {
4713  AssertIndexRange(face, faces_per_cell);
4714  AssertIndexRange(line, lines_per_face);
4715 
4716  const unsigned lines[faces_per_cell][lines_per_face] = {
4717  {8, 10, 0, 4}, // left face
4718  {9, 11, 1, 5}, // right face
4719  {2, 6, 8, 9}, // front face
4720  {3, 7, 10, 11}, // back face
4721  {0, 1, 2, 3}, // bottom face
4722  {4, 5, 6, 7}}; // top face
4723  return lines[face][real_to_standard_face_line(
4724  line, face_orientation, face_flip, face_rotation)];
4725 }
4726 
4727 
4728 
4729 inline unsigned int
4730 GeometryInfo<0>::face_to_cell_lines(const unsigned int,
4731  const unsigned int,
4732  const bool,
4733  const bool,
4734  const bool)
4735 {
4736  Assert(false, ExcNotImplemented());
4738 }
4739 
4740 
4741 
4742 template <int dim>
4743 inline unsigned int
4744 GeometryInfo<dim>::face_to_cell_lines(const unsigned int,
4745  const unsigned int,
4746  const bool,
4747  const bool,
4748  const bool)
4749 {
4750  Assert(false, ExcNotImplemented());
4752 }
4753 
4754 
4755 
4756 template <int dim>
4757 inline unsigned int
4758 GeometryInfo<dim>::face_to_cell_vertices(const unsigned int face,
4759  const unsigned int vertex,
4760  const bool face_orientation,
4761  const bool face_flip,
4762  const bool face_rotation)
4763 {
4764  return child_cell_on_face(RefinementCase<dim>::isotropic_refinement,
4765  face,
4766  vertex,
4767  face_orientation,
4768  face_flip,
4769  face_rotation);
4770 }
4771 
4772 
4773 
4774 inline unsigned int
4775 GeometryInfo<0>::face_to_cell_vertices(const unsigned int,
4776  const unsigned int,
4777  const bool,
4778  const bool,
4779  const bool)
4780 {
4781  Assert(false, ExcNotImplemented());
4783 }
4784 
4785 
4786 
4787 template <int dim>
4788 template <typename Number>
4789 inline Point<dim, Number>
4791 {
4793  for (unsigned int i = 0; i < dim; ++i)
4794  p[i] = std::min(std::max(q[i], Number(0.)), Number(1.));
4795 
4796  return p;
4797 }
4798 
4799 
4800 
4801 template <int dim>
4802 inline double
4804 {
4805  double result = 0.0;
4806 
4807  for (unsigned int i = 0; i < dim; ++i)
4808  {
4809  result = std::max(result, -p[i]);
4810  result = std::max(result, p[i] - 1.);
4811  }
4812 
4813  return result;
4814 }
4815 
4816 
4817 
4818 template <int dim>
4819 inline double
4821  const unsigned int i)
4822 {
4824 
4825  switch (dim)
4826  {
4827  case 1:
4828  {
4829  const double x = xi[0];
4830  switch (i)
4831  {
4832  case 0:
4833  return 1 - x;
4834  case 1:
4835  return x;
4836  }
4837  break;
4838  }
4839 
4840  case 2:
4841  {
4842  const double x = xi[0];
4843  const double y = xi[1];
4844  switch (i)
4845  {
4846  case 0:
4847  return (1 - x) * (1 - y);
4848  case 1:
4849  return x * (1 - y);
4850  case 2:
4851  return (1 - x) * y;
4852  case 3:
4853  return x * y;
4854  }
4855  break;
4856  }
4857 
4858  case 3:
4859  {
4860  const double x = xi[0];
4861  const double y = xi[1];
4862  const double z = xi[2];
4863  switch (i)
4864  {
4865  case 0:
4866  return (1 - x) * (1 - y) * (1 - z);
4867  case 1:
4868  return x * (1 - y) * (1 - z);
4869  case 2:
4870  return (1 - x) * y * (1 - z);
4871  case 3:
4872  return x * y * (1 - z);
4873  case 4:
4874  return (1 - x) * (1 - y) * z;
4875  case 5:
4876  return x * (1 - y) * z;
4877  case 6:
4878  return (1 - x) * y * z;
4879  case 7:
4880  return x * y * z;
4881  }
4882  break;
4883  }
4884 
4885  default:
4886  Assert(false, ExcNotImplemented());
4887  }
4888  return -1e9;
4889 }
4890 
4891 
4892 
4893 template <>
4895  const Point<1> &,
4896  const unsigned int i)
4897 {
4899 
4900  switch (i)
4901  {
4902  case 0:
4903  return Point<1>(-1.);
4904  case 1:
4905  return Point<1>(1.);
4906  }
4907 
4908  return Point<1>(-1e9);
4909 }
4910 
4911 
4912 
4913 template <>
4915  const Point<2> &xi,
4916  const unsigned int i)
4917 {
4919 
4920  const double x = xi[0];
4921  const double y = xi[1];
4922  switch (i)
4923  {
4924  case 0:
4925  return Point<2>(-(1 - y), -(1 - x));
4926  case 1:
4927  return Point<2>(1 - y, -x);
4928  case 2:
4929  return Point<2>(-y, 1 - x);
4930  case 3:
4931  return Point<2>(y, x);
4932  }
4933  return Point<2>(-1e9, -1e9);
4934 }
4935 
4936 
4937 
4938 template <>
4940  const Point<3> &xi,
4941  const unsigned int i)
4942 {
4944 
4945  const double x = xi[0];
4946  const double y = xi[1];
4947  const double z = xi[2];
4948  switch (i)
4949  {
4950  case 0:
4951  return Point<3>(-(1 - y) * (1 - z),
4952  -(1 - x) * (1 - z),
4953  -(1 - x) * (1 - y));
4954  case 1:
4955  return Point<3>((1 - y) * (1 - z), -x * (1 - z), -x * (1 - y));
4956  case 2:
4957  return Point<3>(-y * (1 - z), (1 - x) * (1 - z), -(1 - x) * y);
4958  case 3:
4959  return Point<3>(y * (1 - z), x * (1 - z), -x * y);
4960  case 4:
4961  return Point<3>(-(1 - y) * z, -(1 - x) * z, (1 - x) * (1 - y));
4962  case 5:
4963  return Point<3>((1 - y) * z, -x * z, x * (1 - y));
4964  case 6:
4965  return Point<3>(-y * z, (1 - x) * z, (1 - x) * y);
4966  case 7:
4967  return Point<3>(y * z, x * z, x * y);
4968  }
4969 
4970  return Point<3>(-1e9, -1e9, -1e9);
4971 }
4972 
4973 
4974 
4975 template <int dim>
4976 inline Tensor<1, dim>
4978  const unsigned int)
4979 {
4980  Assert(false, ExcNotImplemented());
4981  return Tensor<1, dim>();
4982 }
4983 
4984 
4985 unsigned int inline GeometryInfo<0>::n_children(const RefinementCase<0> &)
4986 {
4987  return 0;
4988 }
4989 
4990 
4991 namespace internal
4992 {
4993  namespace GeometryInfoHelper
4994  {
4995  // wedge product of a single
4996  // vector in 2d: we just have to
4997  // rotate it by 90 degrees to the
4998  // right
4999  inline Tensor<1, 2>
5000  wedge_product(const Tensor<1, 2> (&derivative)[1])
5001  {
5002  Tensor<1, 2> result;
5003  result[0] = derivative[0][1];
5004  result[1] = -derivative[0][0];
5005 
5006  return result;
5007  }
5008 
5009 
5010  // wedge product of 2 vectors in
5011  // 3d is the cross product
5012  inline Tensor<1, 3>
5013  wedge_product(const Tensor<1, 3> (&derivative)[2])
5014  {
5015  return cross_product_3d(derivative[0], derivative[1]);
5016  }
5017 
5018 
5019  // wedge product of dim vectors
5020  // in dim-d: that's the
5021  // determinant of the matrix
5022  template <int dim>
5023  inline Tensor<0, dim>
5024  wedge_product(const Tensor<1, dim> (&derivative)[dim])
5025  {
5026  Tensor<2, dim> jacobian;
5027  for (unsigned int i = 0; i < dim; ++i)
5028  jacobian[i] = derivative[i];
5029 
5030  return determinant(jacobian);
5031  }
5032  } // namespace GeometryInfoHelper
5033 } // namespace internal
5034 
5035 
5036 template <int dim>
5037 template <int spacedim>
5038 inline void
5040 # ifndef DEAL_II_CXX14_CONSTEXPR_BUG
5041  (const Point<spacedim> (&vertices)[vertices_per_cell],
5042  Tensor<spacedim - dim, spacedim> (&forms)[vertices_per_cell])
5043 # else
5044  (const Point<spacedim> *vertices, Tensor<spacedim - dim, spacedim> *forms)
5045 # endif
5046 {
5047  // for each of the vertices,
5048  // compute the alternating form
5049  // of the mapped unit
5050  // vectors. consider for
5051  // example the case of a quad
5052  // in spacedim==3: to do so, we
5053  // need to see how the
5054  // infinitesimal vectors
5055  // (d\xi_1,0) and (0,d\xi_2)
5056  // are transformed into
5057  // spacedim-dimensional space
5058  // and then form their cross
5059  // product (i.e. the wedge product
5060  // of two vectors). to this end, note
5061  // that
5062  // \vec x = sum_i \vec v_i phi_i(\vec xi)
5063  // so the transformed vectors are
5064  // [x(\xi+(d\xi_1,0))-x(\xi)]/d\xi_1
5065  // and
5066  // [x(\xi+(0,d\xi_2))-x(\xi)]/d\xi_2
5067  // which boils down to the columns
5068  // of the 3x2 matrix \grad_\xi x(\xi)
5069  //
5070  // a similar reasoning would
5071  // hold for all dim,spacedim
5072  // pairs -- we only have to
5073  // compute the wedge product of
5074  // the columns of the
5075  // derivatives
5076  for (unsigned int i = 0; i < vertices_per_cell; ++i)
5077  {
5078  Tensor<1, spacedim> derivatives[dim];
5079 
5080  for (unsigned int j = 0; j < vertices_per_cell; ++j)
5081  {
5082  const Tensor<1, dim> grad_phi_j =
5083  d_linear_shape_function_gradient(unit_cell_vertex(i), j);
5084  for (unsigned int l = 0; l < dim; ++l)
5085  derivatives[l] += vertices[j] * grad_phi_j[l];
5086  }
5087 
5088  forms[i] = internal::GeometryInfoHelper::wedge_product(derivatives);
5089  }
5090 }
5091 
5092 #endif // DOXYGEN
5094 
5095 #endif
GeometryPrimitive(const unsigned int object_dimension)
GeometryPrimitive(const Object object)
Definition: point.h:112
RefinementCase operator|(const RefinementCase &r) const
static std::size_t memory_consumption()
RefinementCase operator~() const
RefinementCase operator&(const RefinementCase &r) const
static std::array< RefinementCase< dim >, n_refinement_cases > all_refinement_cases()
RefinementCase(const typename RefinementPossibilities< dim >::Possibilities refinement_case)
static constexpr unsigned int n_refinement_cases
static RefinementCase cut_axis(const unsigned int i)
std::uint8_t value
void serialize(Archive &ar, const unsigned int version)
RefinementCase(const std::uint8_t refinement_case)
Definition: tensor.h:516
static constexpr std::size_t memory_consumption()
SubfaceCase(const typename SubfacePossibilities< dim >::Possibilities subface_possibility)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
Point< 3 > vertices[4]
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcInvalidCoordinate(double arg1)
static ::ExceptionBase & ExcInvalidSubface(int arg1, int arg2, int arg3)
#define Assert(cond, exc)
Definition: exceptions.h:1631
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcInvalidSubfaceCase(int arg1)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1888
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition: exceptions.h:563
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:517
static ::ExceptionBase & ExcInvalidRefinementCase(int arg1)
static ::ExceptionBase & ExcMessage(std::string arg1)
static const char U
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:192
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static const unsigned int invalid_unsigned_int
Definition: types.h:221
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
typename internal::ndarray::HelperArray< T, Ns... >::type ndarray
Definition: ndarray.h:108
#define DEAL_II_HOST_DEVICE
Definition: numbers.h:35
static unsigned int n_children(const RefinementCase< 0 > &refinement_case)
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 const std::array< unsigned int, vertices_per_cell > dx_to_deal
static std::array< unsigned int, vertices_per_cell > vertex_indices()
static std::array< unsigned int, 0 > face_indices()
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 const std::array< unsigned int, vertices_per_cell > ucd_to_deal
static unsigned int child_cell_from_point(const Point< dim > &p)
static unsigned int standard_to_real_line_vertex(const unsigned int vertex, const bool line_orientation=true)
static constexpr std::array< unsigned int, vertices_per_cell > ucd_to_deal
static constexpr ndarray< Tensor< 1, dim >, faces_per_cell, dim - 1 > unit_tangential_vectors
static bool is_inside_unit_cell(const Point< dim > &p, const double eps)
static RefinementCase< dim > min_cell_refinement_case_for_face_refinement(const RefinementCase< dim - 1 > &face_refinement_case, const unsigned int face_no, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
static double distance_to_unit_cell(const Point< dim > &p)
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 constexpr unsigned int quads_per_face
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)
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices()
static constexpr std::array< unsigned int, faces_per_cell > opposite_face
static constexpr unsigned int max_children_per_face
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 constexpr unsigned int vertices_per_cell
static std::array< unsigned int, 2 > standard_hex_line_to_quad_line_index(const unsigned int line)
static Point< dim > unit_cell_vertex(const unsigned int vertex)
static constexpr std::array< unsigned int, vertices_per_cell > dx_to_deal
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices()
static constexpr unsigned int lines_per_cell
static RefinementCase< 1 > line_refinement_case(const RefinementCase< dim > &cell_refinement_case, const unsigned int line_no)
static constexpr std::array< unsigned int, faces_per_cell > unit_normal_direction
static RefinementCase< dim > min_cell_refinement_case_for_line_refinement(const unsigned int line_no)
static constexpr unsigned int faces_per_cell
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 subface_ratio(const internal::SubfaceCase< dim > &subface_case, const unsigned int subface_no)
static constexpr unsigned int hexes_per_cell
static Point< dim > child_to_cell_coordinates(const Point< dim > &p, const unsigned int child_index, const RefinementCase< dim > refine_case=RefinementCase< dim >::isotropic_refinement)
static Point< dim > cell_to_child_coordinates(const Point< dim > &p, const unsigned int child_index, const RefinementCase< dim > refine_case=RefinementCase< dim >::isotropic_refinement)
static double d_linear_shape_function(const Point< dim > &xi, const unsigned int i)
static constexpr std::array< Tensor< 1, dim >, faces_per_cell > unit_normal_vector
static Point< dim, Number > project_to_unit_cell(const Point< dim, Number > &p)
static Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i)
static constexpr unsigned int vertices_per_face
static unsigned int real_to_standard_face_line(const unsigned int line, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
static unsigned int line_to_cell_vertices(const unsigned int line, const unsigned int vertex)
static bool is_inside_unit_cell(const Point< dim > &p)
static unsigned int n_children(const RefinementCase< dim > &refinement_case)
static unsigned int real_to_standard_face_vertex(const unsigned int vertex, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
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 constexpr ndarray< unsigned int, vertices_per_cell, dim > vertex_to_face
static void alternating_form_at_vertices(const Point< spacedim >(&vertices)[vertices_per_cell], Tensor< spacedim - dim, spacedim >(&forms)[vertices_per_cell])
static RefinementCase< dim - 1 > face_refinement_case(const RefinementCase< dim > &cell_refinement_case, const unsigned int face_no, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
static constexpr std::array< int, faces_per_cell > unit_normal_orientation
static constexpr unsigned int quads_per_cell
static constexpr unsigned int max_children_per_cell
static unsigned int n_subfaces(const internal::SubfaceCase< dim > &subface_case)
static constexpr unsigned int lines_per_face
static std::array< unsigned int, 2 > standard_hex_vertex_to_quad_vertex_index(const unsigned int vertex)
constexpr DEAL_II_HOST Number determinant(const SymmetricTensor< 2, dim, Number > &)