Reference documentation for deal.II version GIT f6a5d312c9 2023-10-04 08:50:02+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\}}\)
tria_accessor.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_tria_accessor_h
17 #define dealii_tria_accessor_h
18 
19 
20 #include <deal.II/base/config.h>
21 
25 #include <deal.II/base/point.h>
26 
27 #include <deal.II/grid/cell_id.h>
31 
32 #include <boost/container/small_vector.hpp>
33 
34 #include <utility>
35 
36 
38 
39 // Forward declarations
40 #ifndef DOXYGEN
41 template <int dim, int spacedim>
42 DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
43 class Triangulation;
44 template <typename Accessor>
45 class TriaRawIterator;
46 template <typename Accessor>
47 class TriaIterator;
48 template <typename Accessor>
49 class TriaActiveIterator;
50 
51 namespace parallel
52 {
53  template <int dim, int spacedim>
54  DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
55  class TriangulationBase;
56 }
57 
58 template <int dim, int spacedim>
59 DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
60 class DoFHandler;
61 template <int dim, int spacedim, bool lda>
62 class DoFCellAccessor;
63 
64 
65 template <int dim, int spacedim>
66 class Manifold;
67 
68 template <int dim, int spacedim>
69 class Mapping;
70 #endif
71 
72 namespace internal
73 {
74  namespace TriangulationImplementation
75  {
76  class TriaObjects;
77  struct Implementation;
78  struct ImplementationMixedMesh;
79  } // namespace TriangulationImplementation
80 
81  namespace TriaAccessorImplementation
82  {
83  struct Implementation;
84 
90  template <int structdim, int dim>
92  {
93  struct type
94  {
98  type() = default;
99 
103  type(const int level)
104  {
105  Assert(level == 0, ExcInternalError());
106  (void)level; // removes -Wunused-parameter warning in optimized mode
107  }
108 
112  operator int() const
113  {
114  return 0;
115  }
116 
117  void
118  operator++() const
119  {
120  Assert(false, ExcInternalError());
121  }
122 
123  void
124  operator--() const
125  {
126  Assert(false, ExcInternalError());
127  }
128  };
129  };
130 
131 
137  template <int dim>
138  struct PresentLevelType<dim, dim>
139  {
140  using type = int;
141  };
142  } // namespace TriaAccessorImplementation
143 } // namespace internal
144 template <int structdim, int dim, int spacedim>
145 class TriaAccessor;
146 template <int dim, int spacedim>
147 class TriaAccessor<0, dim, spacedim>;
148 template <int spacedim>
149 class TriaAccessor<0, 1, spacedim>;
150 
151 // note: the file tria_accessor.templates.h is included at the end of
152 // this file. this includes a lot of templates. originally, this was
153 // only done in debug mode, but led to cyclic reduction problems and
154 // so is now on by default.
155 
156 
161 {
166  "The operation you are attempting can only be performed for "
167  "(cell, face, or edge) iterators that point to valid "
168  "objects. These objects need not necessarily be active, "
169  "i.e., have no children, but they need to be part of a "
170  "triangulation. (The objects pointed to by an iterator "
171  "may -- after coarsening -- also be objects that used "
172  "to be part of a triangulation, but are now no longer "
173  "used. Their memory location may have been retained "
174  "for re-use upon the next mesh refinement, but is "
175  "currently unused.)");
186  "The operation you are attempting can only be performed for "
187  "(cell, face, or edge) iterators that point to 'active' "
188  "objects. 'Active' objects are those that do not have "
189  "children (in the case of cells), or that are part of "
190  "an active cell (in the case of faces or edges). However, "
191  "the object on which you are trying the current "
192  "operation is not 'active' in this sense.");
199  "The operation you are attempting can only be performed for "
200  "(cell, face, or edge) iterators that have children, "
201  "but the object on which you are trying the current "
202  "operation does not have any.");
210  "The operation you are attempting can only be performed for "
211  "(cell, face, or edge) iterators that have a parent object, "
212  "but the object on which you are trying the current "
213  "operation does not have one -- i.e., it is on the "
214  "coarsest level of the triangulation.");
219  int,
220  << "You can only set the child index if the cell does not "
221  << "currently have children registered; or you can clear it. "
222  << "The given index was " << arg1
223  << " (-1 means: clear children).");
227  template <typename AccessorType>
229  AccessorType,
230  << "You tried to dereference an iterator for which this "
231  << "is not possible. More information on this iterator: "
232  << "index=" << arg1.index() << ", state="
233  << (arg1.state() == IteratorState::valid ?
234  "valid" :
235  (arg1.state() == IteratorState::past_the_end ?
236  "past_the_end" :
237  "invalid")));
242  "Iterators can only be compared if they point to the same "
243  "triangulation, or if neither of them are associated "
244  "with a triangulation.");
245  // TODO: Write documentation!
250  // TODO: Write documentation!
270  // TODO: Write documentation!
276  int,
277  << "You can only set the child index of an even numbered child."
278  << "The number of the child given was " << arg1 << '.');
279 } // namespace TriaAccessorExceptions
280 
281 
307 template <int structdim, int dim, int spacedim = dim>
309 {
310 public:
317  static constexpr unsigned int space_dimension = spacedim;
318 
324  static constexpr unsigned int dimension = dim;
325 
331  static const unsigned int structure_dimension = structdim;
332 
342  void
343  operator=(const TriaAccessorBase *) = delete;
344 
345 protected:
351  using AccessorData = void;
352 
357  const int level = -1,
358  const int index = -1,
359  const AccessorData * = nullptr);
360 
365 
373  void
375 
381 
392  bool
393  operator<(const TriaAccessorBase &other) const;
394 
395 protected:
399  bool
400  operator==(const TriaAccessorBase &) const;
401 
405  bool
406  operator!=(const TriaAccessorBase &) const;
407 
421  void
423 
431  void
441  objects() const;
442 
443 public:
449  using LocalData = void *;
450 
474  int
475  level() const;
476 
503  int
504  index() const;
505 
511  state() const;
512 
519 
523 protected:
528  typename ::internal::TriaAccessorImplementation::
529  PresentLevelType<structdim, dim>::type present_level;
530 
536 
541 
542 private:
543  template <typename Accessor>
544  friend class TriaRawIterator;
545  template <typename Accessor>
546  friend class TriaIterator;
547  template <typename Accessor>
548  friend class TriaActiveIterator;
549 };
550 
551 
552 
573 template <int structdim, int dim, int spacedim = dim>
575 {
576 public:
583  static constexpr unsigned int space_dimension = spacedim;
584 
590  static constexpr unsigned int dimension = dim;
591 
597  static const unsigned int structure_dimension = structdim;
598 
604  using AccessorData = void;
605 
613  InvalidAccessor(const void *parent = nullptr,
614  const int level = -1,
615  const int index = -1,
616  const AccessorData *local_data = nullptr);
617 
626 
631  template <typename OtherAccessor>
632  InvalidAccessor(const OtherAccessor &);
633 
637  void
639 
643  bool
644  operator==(const InvalidAccessor &) const;
645  bool
646  operator!=(const InvalidAccessor &) const;
647 
651  void
652  operator++() const;
653  void
654  operator--() const;
655 
661  state();
662 
663 
668  static int
669  level();
670 
675  static int
676  index();
677 
682  bool
683  used() const;
684 
689  bool
690  has_children() const;
691 
696  manifold_id() const;
697 
701  unsigned int
702  user_index() const;
703 
707  void
708  set_user_index(const unsigned int p) const;
709 
713  void
715 
720  vertex(const unsigned int i) const;
721 
726  void *
727  line(const unsigned int i) const;
728 
733  void *
734  quad(const unsigned int i) const;
735 };
736 
737 
738 
756 template <int structdim, int dim, int spacedim>
757 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
758 {
759 public:
763  using AccessorData =
765 
770  const int level = -1,
771  const int index = -1,
772  const AccessorData *local_data = nullptr);
773 
778  TriaAccessor(const TriaAccessor &) = default;
779 
783  TriaAccessor(TriaAccessor &&) = default; // NOLINT
784 
797  template <int structdim2, int dim2, int spacedim2>
799 
804  template <int structdim2, int dim2, int spacedim2>
806 
816  TriaAccessor &
817  operator=(const TriaAccessor &) = delete;
818 
822  TriaAccessor &
823  operator=(TriaAccessor &&) = default; // NOLINT
824 
828  ~TriaAccessor() = default;
829 
836  bool
837  used() const;
838 
851  vertex_iterator(const unsigned int i) const;
852 
868  unsigned int
869  vertex_index(const unsigned int i) const;
870 
909  vertex(const unsigned int i) const;
910 
914  typename ::internal::TriangulationImplementation::
915  Iterators<dim, spacedim>::line_iterator
916  line(const unsigned int i) const;
917 
924  unsigned int
925  line_index(const unsigned int i) const;
926 
930  typename ::internal::TriangulationImplementation::
931  Iterators<dim, spacedim>::quad_iterator
932  quad(const unsigned int i) const;
933 
940  unsigned int
941  quad_index(const unsigned int i) const;
959  unsigned char
960  combined_face_orientation(const unsigned int face) const;
961 
973  bool
974  face_orientation(const unsigned int face) const;
975 
985  bool
986  face_flip(const unsigned int face) const;
987 
997  bool
998  face_rotation(const unsigned int face) const;
999 
1012  bool
1013  line_orientation(const unsigned int line) const;
1028  bool
1029  has_children() const;
1030 
1035  unsigned int
1036  n_children() const;
1037 
1051  unsigned int
1053 
1067  unsigned int
1069 
1074  child(const unsigned int i) const;
1075 
1080  unsigned int
1083 
1093  isotropic_child(const unsigned int i) const;
1094 
1100 
1106  int
1107  child_index(const unsigned int i) const;
1108 
1114  int
1115  isotropic_child_index(const unsigned int i) const;
1138  boundary_id() const;
1139 
1169  void
1171 
1202  void
1204 
1212  bool
1213  at_boundary() const;
1214 
1224  const Manifold<dim, spacedim> &
1225  get_manifold() const;
1226 
1248  manifold_id() const;
1249 
1267  void
1269 
1283  void
1285 
1302  bool
1303  user_flag_set() const;
1304 
1310  void
1311  set_user_flag() const;
1312 
1318  void
1320 
1326  void
1328 
1334  void
1336 
1342  void
1344 
1356  void
1357  set_user_pointer(void *p) const;
1358 
1364  void
1366 
1382  void *
1383  user_pointer() const;
1384 
1406  void
1408 
1415  void
1417 
1427  void
1428  set_user_index(const unsigned int p) const;
1429 
1435  void
1437 
1449  unsigned int
1450  user_index() const;
1451 
1469  void
1470  recursively_set_user_index(const unsigned int p) const;
1471 
1480  void
1516  double
1517  diameter() const;
1518 
1545  std::pair<Point<spacedim>, double>
1547 
1557  bounding_box() const;
1558 
1568  double
1569  extent_in_direction(const unsigned int axis) const;
1570 
1574  double
1576 
1591  intermediate_point(const Point<structdim> &coordinates) const;
1592 
1617 
1653  center(const bool respect_manifold = false,
1654  const bool interpolate_from_surrounding = false) const;
1655 
1674  barycenter() const;
1675 
1706  double
1707  measure() const;
1708 
1723  bool
1726 
1732 
1736  unsigned int
1737  n_vertices() const;
1738 
1742  unsigned int
1743  n_lines() const;
1744 
1754  unsigned int
1755  n_faces() const;
1756 
1763 
1769  line_indices() const;
1770 
1778  face_indices() const;
1779 
1784 private:
1789  void
1791 
1799  void
1801  const std::initializer_list<int> &new_indices) const;
1802 
1806  void
1808  const std::initializer_list<unsigned int> &new_indices) const;
1809 
1817  void
1818  set_line_orientation(const unsigned int line, const bool orientation) const;
1819 
1829  void
1830  set_combined_face_orientation(const unsigned int face,
1831  const unsigned char combined_orientation) const;
1832 
1836  void
1837  set_used_flag() const;
1838 
1842  void
1844 
1853  void
1855 
1863  void
1865 
1872  void
1873  set_children(const unsigned int i, const int index) const;
1874 
1879  void
1881 
1882 private:
1883  friend class Triangulation<dim, spacedim>;
1884 
1885  friend struct ::internal::TriangulationImplementation::Implementation;
1886  friend struct ::internal::TriangulationImplementation::
1887  ImplementationMixedMesh;
1888  friend struct ::internal::TriaAccessorImplementation::Implementation;
1889 };
1890 
1891 
1892 
1911 template <int dim, int spacedim>
1912 class TriaAccessor<0, dim, spacedim>
1913 {
1914 public:
1920  static constexpr unsigned int space_dimension = spacedim;
1921 
1927  static constexpr unsigned int dimension = dim;
1928 
1934  static const unsigned int structure_dimension = 0;
1935 
1939  using AccessorData = void;
1940 
1946  const unsigned int vertex_index);
1947 
1954  const int level = 0,
1955  const int index = 0,
1956  const AccessorData * = nullptr);
1957 
1961  template <int structdim2, int dim2, int spacedim2>
1963 
1967  template <int structdim2, int dim2, int spacedim2>
1969 
1974  state() const;
1975 
1980  static int
1982 
1987  int
1988  index() const;
1989 
1996 
2006  void
2008 
2012  void
2017  bool
2018  operator==(const TriaAccessor &) const;
2019 
2023  bool
2024  operator!=(const TriaAccessor &) const;
2025 
2053  unsigned int
2054  vertex_index(const unsigned int i = 0) const;
2055 
2061  Point<spacedim> &
2062  vertex(const unsigned int i = 0) const;
2063 
2068  typename ::internal::TriangulationImplementation::
2069  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2070 
2074  static unsigned int
2075  line_index(const unsigned int i);
2076 
2080  static typename ::internal::TriangulationImplementation::
2081  Iterators<dim, spacedim>::quad_iterator
2082  quad(const unsigned int i);
2083 
2087  static unsigned int
2088  quad_index(const unsigned int i);
2089 
2105  double
2106  diameter() const;
2107 
2115  double
2116  extent_in_direction(const unsigned int axis) const;
2117 
2126  center(const bool respect_manifold = false,
2127  const bool interpolate_from_surrounding = false) const;
2128 
2137  double
2138  measure() const;
2155  static unsigned char
2156  combined_face_orientation(const unsigned int face);
2157 
2161  static bool
2162  face_orientation(const unsigned int face);
2163 
2167  static bool
2168  face_flip(const unsigned int face);
2169 
2173  static bool
2174  face_rotation(const unsigned int face);
2175 
2179  static bool
2180  line_orientation(const unsigned int line);
2181 
2196  static bool
2198 
2203  static unsigned int
2205 
2210  static unsigned int
2212 
2216  static unsigned int
2218 
2222  static unsigned int
2224 
2229  child(const unsigned int);
2230 
2235  isotropic_child(const unsigned int);
2236 
2240  static RefinementCase<0>
2242 
2246  static int
2247  child_index(const unsigned int i);
2248 
2252  static int
2253  isotropic_child_index(const unsigned int i);
2261  bool
2262  used() const;
2263 
2264 protected:
2272  void
2274 
2283  bool
2284  operator<(const TriaAccessor &other) const;
2285 
2290 
2294  unsigned int global_vertex_index;
2295 
2296 private:
2297  template <typename Accessor>
2298  friend class TriaRawIterator;
2299  template <typename Accessor>
2300  friend class TriaIterator;
2301  template <typename Accessor>
2302  friend class TriaActiveIterator;
2303 };
2304 
2305 
2306 
2323 template <int spacedim>
2324 class TriaAccessor<0, 1, spacedim>
2325 {
2326 public:
2332  static constexpr unsigned int space_dimension = spacedim;
2333 
2339  static constexpr unsigned int dimension = 1;
2340 
2346  static const unsigned int structure_dimension = 0;
2347 
2351  using AccessorData = void;
2352 
2358  {
2370  right_vertex
2371  };
2372 
2385  const VertexKind vertex_kind,
2386  const unsigned int vertex_index);
2387 
2394  const int = 0,
2395  const int = 0,
2396  const AccessorData * = nullptr);
2397 
2401  template <int structdim2, int dim2, int spacedim2>
2403 
2407  template <int structdim2, int dim2, int spacedim2>
2409 
2414  void
2416 
2422  void
2424 
2432 
2437  static int
2439 
2444  int
2445  index() const;
2446 
2453 
2464  void
2465  operator++() const;
2466 
2471  void
2472  operator--() const;
2476  bool
2477  operator==(const TriaAccessor &) const;
2478 
2482  bool
2483  operator!=(const TriaAccessor &) const;
2484 
2493  bool
2494  operator<(const TriaAccessor &other) const;
2495 
2522  unsigned int
2523  vertex_index(const unsigned int i = 0) const;
2524 
2530  Point<spacedim> &
2531  vertex(const unsigned int i = 0) const;
2532 
2538  center() const;
2539 
2544  typename ::internal::TriangulationImplementation::
2545  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2546 
2553  static unsigned int
2554  line_index(const unsigned int i);
2555 
2559  static typename ::internal::TriangulationImplementation::
2560  Iterators<1, spacedim>::quad_iterator
2561  quad(const unsigned int i);
2562 
2569  static unsigned int
2570  quad_index(const unsigned int i);
2571 
2581  bool
2582  at_boundary() const;
2583 
2599  boundary_id() const;
2600 
2604  const Manifold<1, spacedim> &
2605  get_manifold() const;
2606 
2614  manifold_id() const;
2615 
2616 
2628  bool
2629  user_flag_set() const;
2630 
2636  void
2637  set_user_flag() const;
2638 
2644  void
2645  clear_user_flag() const;
2646 
2652  void
2653  recursively_set_user_flag() const;
2654 
2660  void
2662 
2668  void
2669  clear_user_data() const;
2670 
2682  void
2683  set_user_pointer(void *p) const;
2684 
2690  void
2691  clear_user_pointer() const;
2692 
2708  void *
2709  user_pointer() const;
2710 
2732  void
2733  recursively_set_user_pointer(void *p) const;
2734 
2741  void
2743 
2753  void
2754  set_user_index(const unsigned int p) const;
2755 
2761  void
2762  clear_user_index() const;
2763 
2775  unsigned int
2776  user_index() const;
2777 
2795  void
2796  recursively_set_user_index(const unsigned int p) const;
2797 
2806  void
2824  static unsigned char
2825  combined_face_orientation(const unsigned int face);
2826 
2830  static bool
2831  face_orientation(const unsigned int face);
2832 
2836  static bool
2837  face_flip(const unsigned int face);
2838 
2842  static bool
2843  face_rotation(const unsigned int face);
2844 
2848  static bool
2849  line_orientation(const unsigned int line);
2850 
2865  static bool
2867 
2872  static unsigned int
2874 
2879  static unsigned int
2881 
2885  static unsigned int
2887 
2891  static unsigned int
2893 
2898  child(const unsigned int);
2899 
2904  isotropic_child(const unsigned int);
2905 
2909  static RefinementCase<0>
2911 
2915  static int
2916  child_index(const unsigned int i);
2917 
2921  static int
2922  isotropic_child_index(const unsigned int i);
2955  void
2957 
2964  void
2966 
2978  void
2980 
2992  void
3001  bool
3002  used() const;
3003 
3009 
3013  unsigned int
3014  n_vertices() const;
3015 
3019  unsigned int
3020  n_lines() const;
3021 
3028 
3034  line_indices() const;
3035 
3036 protected:
3041 
3047 
3051  unsigned int global_vertex_index;
3052 };
3053 
3054 
3055 
3071 template <int dim, int spacedim = dim>
3072 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
3073 {
3074 public:
3079 
3084 
3096  const int level = -1,
3097  const int index = -1,
3098  const AccessorData *local_data = nullptr);
3099 
3104 
3117  template <int structdim2, int dim2, int spacedim2>
3119 
3124  template <int structdim2, int dim2, int spacedim2>
3126 
3131 
3135  // NOLINTNEXTLINE OSX does not compile with noexcept
3137 
3141  ~CellAccessor() = default;
3142 
3154 
3158  // NOLINTNEXTLINE OSX does not compile with noexcept
3160  operator=(CellAccessor<dim, spacedim> &&) = default; // NOLINT
3161 
3185  as_dof_handler_iterator(const DoFHandler<dim, spacedim> &dof_handler) const;
3186 
3197  const DoFHandler<dim, spacedim> &dof_handler) const;
3198 
3199 
3216  child(const unsigned int i) const;
3217 
3221  boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
3224 
3228  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
3229  face(const unsigned int i) const;
3230 
3235  unsigned int
3238 
3242  boost::container::small_vector<
3243  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
3246 
3256  unsigned int
3257  face_index(const unsigned int i) const;
3258 
3307  neighbor_child_on_subface(const unsigned int face_no,
3308  const unsigned int subface_no) const;
3309 
3335  neighbor(const unsigned int face_no) const;
3336 
3344  int
3345  neighbor_index(const unsigned int face_no) const;
3346 
3354  int
3355  neighbor_level(const unsigned int face_no) const;
3356 
3368  unsigned int
3369  neighbor_of_neighbor(const unsigned int face_no) const;
3370 
3381  bool
3382  neighbor_is_coarser(const unsigned int face_no) const;
3383 
3398  std::pair<unsigned int, unsigned int>
3399  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3400 
3407  unsigned int
3408  neighbor_face_no(const unsigned int neighbor) const;
3409 
3413  static bool
3415 
3429  bool
3430  has_periodic_neighbor(const unsigned int i) const;
3431 
3449  periodic_neighbor(const unsigned int i) const;
3450 
3459  neighbor_or_periodic_neighbor(const unsigned int i) const;
3460 
3476  periodic_neighbor_child_on_subface(const unsigned int face_no,
3477  const unsigned int subface_no) const;
3478 
3489  std::pair<unsigned int, unsigned int>
3490  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3491 
3497  int
3498  periodic_neighbor_index(const unsigned int i) const;
3499 
3505  int
3506  periodic_neighbor_level(const unsigned int i) const;
3507 
3522  unsigned int
3523  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3524 
3530  unsigned int
3531  periodic_neighbor_face_no(const unsigned int i) const;
3532 
3539  bool
3540  periodic_neighbor_is_coarser(const unsigned int i) const;
3541 
3558  bool
3559  at_boundary(const unsigned int i) const;
3560 
3569  bool
3570  at_boundary() const;
3571 
3579  bool
3580  has_boundary_lines() const;
3608 
3626  void
3629 
3633  void
3635 
3643  bool
3645  const unsigned int face_no,
3646  const RefinementCase<dim - 1> &face_refinement_case =
3648 
3654  bool
3655  flag_for_line_refinement(const unsigned int line_no) const;
3656 
3666  subface_case(const unsigned int face_no) const;
3667 
3671  bool
3673 
3678  void
3680 
3684  void
3709  material_id() const;
3710 
3722  void
3723  set_material_id(const types::material_id new_material_id) const;
3724 
3733  void
3734  recursively_set_material_id(const types::material_id new_material_id) const;
3761  subdomain_id() const;
3762 
3778  void
3779  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3780 
3791 
3796  void
3798  const types::subdomain_id new_level_subdomain_id) const;
3799 
3800 
3816  void
3818  const types::subdomain_id new_subdomain_id) const;
3842 
3852 
3866  bool
3867  direction_flag() const;
3868 
3894  unsigned int
3896 
3904  int
3905  parent_index() const;
3906 
3913  parent() const;
3914 
3934  bool
3935  is_active() const;
3936 
3956  bool
3958 
3963  bool
3965 
3999  bool
4000  is_ghost() const;
4001 
4007  bool
4009 
4036  bool
4037  is_artificial() const;
4038 
4045  bool
4047 
4061  bool
4063 
4072  void
4073  set_neighbor(const unsigned int i,
4074  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
4075 
4089  CellId
4090  id() const;
4091 
4093 
4097  double
4098  diameter(const Mapping<dim, spacedim> &mapping) const;
4099 
4117 
4118 protected:
4134  unsigned int
4135  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
4136 
4142  template <int dim_, int spacedim_>
4143  bool
4144  point_inside_codim(const Point<spacedim_> &p) const;
4145 
4146 
4147 
4148 private:
4153  void
4154  set_active_cell_index(const unsigned int active_cell_index) const;
4155 
4159  void
4161 
4165  void
4167 
4171  void
4172  set_parent(const unsigned int parent_index);
4173 
4180  void
4181  set_direction_flag(const bool new_direction_flag) const;
4182 
4183  friend class Triangulation<dim, spacedim>;
4184 
4185  friend class parallel::TriangulationBase<dim, spacedim>;
4186 
4187  friend struct ::internal::TriangulationImplementation::Implementation;
4188  friend struct ::internal::TriangulationImplementation::
4189  ImplementationMixedMesh;
4190 };
4191 
4192 
4193 
4194 /* ----- declaration of explicit specializations and general templates ----- */
4195 
4196 
4197 template <int structdim, int dim, int spacedim>
4198 template <typename OtherAccessor>
4200  const OtherAccessor &)
4201 {
4202  Assert(false,
4203  ExcMessage("You are attempting an illegal conversion between "
4204  "iterator/accessor types. The constructor you call "
4205  "only exists to make certain template constructs "
4206  "easier to write as dimension independent code but "
4207  "the conversion is not valid in the current context."));
4208 }
4209 
4210 
4211 
4212 template <int structdim, int dim, int spacedim>
4213 template <int structdim2, int dim2, int spacedim2>
4216 {
4217  Assert(false,
4218  ExcMessage("You are attempting an illegal conversion between "
4219  "iterator/accessor types. The constructor you call "
4220  "only exists to make certain template constructs "
4221  "easier to write as dimension independent code but "
4222  "the conversion is not valid in the current context."));
4223 }
4224 
4225 
4226 
4227 template <int dim, int spacedim>
4228 template <int structdim2, int dim2, int spacedim2>
4231 {
4232  Assert(false,
4233  ExcMessage("You are attempting an illegal conversion between "
4234  "iterator/accessor types. The constructor you call "
4235  "only exists to make certain template constructs "
4236  "easier to write as dimension independent code but "
4237  "the conversion is not valid in the current context."));
4238 }
4239 
4240 
4241 
4242 template <int structdim, int dim, int spacedim>
4243 template <int structdim2, int dim2, int spacedim2>
4246 {
4247  Assert(false,
4248  ExcMessage("You are attempting an illegal conversion between "
4249  "iterator/accessor types. The constructor you call "
4250  "only exists to make certain template constructs "
4251  "easier to write as dimension independent code but "
4252  "the conversion is not valid in the current context."));
4253 }
4254 
4255 
4256 
4257 template <int dim, int spacedim>
4258 template <int structdim2, int dim2, int spacedim2>
4261 {
4262  Assert(false,
4263  ExcMessage("You are attempting an illegal conversion between "
4264  "iterator/accessor types. The constructor you call "
4265  "only exists to make certain template constructs "
4266  "easier to write as dimension independent code but "
4267  "the conversion is not valid in the current context."));
4268 }
4269 
4270 
4271 #ifndef DOXYGEN
4272 
4273 template <>
4274 bool
4276 template <>
4277 bool
4279 template <>
4280 bool
4282 template <>
4283 bool
4285 template <>
4286 bool
4288 template <>
4289 bool
4291 // -------------------------------------------------------------------
4292 
4293 template <>
4294 void
4296 
4297 #endif // DOXYGEN
4298 
4300 
4301 // include more templates in debug and optimized mode
4302 #include <deal.II/grid/tria_accessor.templates.h>
4303 
4304 #endif
void recursively_set_subdomain_id(const types::subdomain_id new_subdomain_id) const
CellAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
TriaIterator< CellAccessor< dim, spacedim > > parent() const
unsigned int neighbor_face_no(const unsigned int neighbor) const
unsigned int face_iterator_to_index(const TriaIterator< TriaAccessor< dim - 1, dim, spacedim >> &face) const
void set_active_cell_index(const unsigned int active_cell_index) const
unsigned int neighbor_of_neighbor_internal(const unsigned int neighbor) const
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor(const unsigned int i) const
types::global_cell_index global_active_cell_index() const
bool is_artificial_on_level() const
void set_direction_flag(const bool new_direction_flag) const
void recursively_set_material_id(const types::material_id new_material_id) const
void set_level_subdomain_id(const types::subdomain_id new_level_subdomain_id) const
types::subdomain_id level_subdomain_id() const
boost::container::small_vector< TriaIterator< CellAccessor< dim, spacedim > >, GeometryInfo< dim >::max_children_per_cell > child_iterators() const
bool flag_for_face_refinement(const unsigned int face_no, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement) const
unsigned int face_index(const unsigned int i) const
double diameter(const Mapping< dim, spacedim > &mapping) const
TriaActiveIterator< DoFCellAccessor< dim, spacedim, false > > as_dof_handler_iterator(const DoFHandler< dim, spacedim > &dof_handler) const
::internal::SubfaceCase< dim > subface_case(const unsigned int face_no) const
bool is_active() const
bool is_ghost() const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
void set_subdomain_id(const types::subdomain_id new_subdomain_id) const
CellAccessor< dim, spacedim > & operator=(CellAccessor< dim, spacedim > &&)=default
bool neighbor_is_coarser(const unsigned int face_no) const
void set_global_level_cell_index(const types::global_cell_index index) const
bool has_periodic_neighbor(const unsigned int i) const
int periodic_neighbor_level(const unsigned int i) const
std::pair< unsigned int, unsigned int > neighbor_of_coarser_neighbor(const unsigned int neighbor) const
~CellAccessor()=default
CellAccessor(const CellAccessor< dim, spacedim > &)=default
void set_coarsen_flag() const
TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > face(const unsigned int i) const
unsigned int neighbor_of_neighbor(const unsigned int face_no) const
void set_material_id(const types::material_id new_material_id) const
bool is_locally_owned() const
TriaIterator< CellAccessor< dim, spacedim > > neighbor(const unsigned int face_no) const
void set_refine_flag(const RefinementCase< dim > ref_case=RefinementCase< dim >::isotropic_refinement) const
RefinementCase< dim > refine_flag_set() const
CellAccessor(CellAccessor< dim, spacedim > &&)=default
bool point_inside_codim(const Point< spacedim_ > &p) const
typename TriaAccessor< dim, dim, spacedim >::AccessorData AccessorData
bool is_locally_owned_on_level() const
bool has_boundary_lines() const
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
int neighbor_level(const unsigned int face_no) const
int periodic_neighbor_index(const unsigned int i) const
bool periodic_neighbor_is_coarser(const unsigned int i) const
void set_global_active_cell_index(const types::global_cell_index index) const
void clear_coarsen_flag() const
TriaIterator< CellAccessor< dim, spacedim > > child(const unsigned int i) const
void set_neighbor(const unsigned int i, const TriaIterator< CellAccessor< dim, spacedim >> &pointer) const
void set_parent(const unsigned int parent_index)
std::pair< unsigned int, unsigned int > periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const
CellAccessor(const TriaAccessor< dim, dim, spacedim > &cell_accessor)
bool at_boundary() const
unsigned int active_cell_index() const
static bool is_level_cell()
void clear_refine_flag() const
int neighbor_index(const unsigned int face_no) const
bool point_inside(const Point< spacedim > &p) const
types::subdomain_id subdomain_id() const
bool direction_flag() const
types::material_id material_id() const
bool coarsen_flag_set() const
types::global_cell_index global_level_cell_index() const
bool flag_for_line_refinement(const unsigned int line_no) const
boost::container::small_vector< TriaIterator< TriaAccessor< dim - 1, dim, spacedim > >, GeometryInfo< dim >::faces_per_cell > face_iterators() const
CellId id() const
bool is_artificial() const
TriaIterator< DoFCellAccessor< dim, spacedim, true > > as_dof_handler_level_iterator(const DoFHandler< dim, spacedim > &dof_handler) const
bool is_ghost_on_level() const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_or_periodic_neighbor(const unsigned int i) const
int parent_index() const
unsigned int periodic_neighbor_of_periodic_neighbor(const unsigned int i) const
unsigned int periodic_neighbor_face_no(const unsigned int i) const
CellAccessor< dim, spacedim > & operator=(const CellAccessor< dim, spacedim > &)=delete
Definition: cell_id.h:72
void * quad(const unsigned int i) const
void set_manifold_id(const types::manifold_id) const
InvalidAccessor(const InvalidAccessor &)
static constexpr unsigned int space_dimension
static int level()
bool operator==(const InvalidAccessor &) const
void operator++() const
static IteratorState::IteratorStates state()
bool operator!=(const InvalidAccessor &) const
bool has_children() const
unsigned int user_index() const
Point< spacedim > & vertex(const unsigned int i) const
static const unsigned int structure_dimension
void copy_from(const InvalidAccessor &)
void set_user_index(const unsigned int p) const
static int index()
types::manifold_id manifold_id() const
InvalidAccessor(const void *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
static constexpr unsigned int dimension
void operator--() const
bool used() const
InvalidAccessor(const OtherAccessor &)
void * line(const unsigned int i) const
Abstract base class for mapping classes.
Definition: mapping.h:317
bool operator!=(const TriaAccessorBase &) const
static constexpr unsigned int space_dimension
static constexpr unsigned int dimension
TriaAccessorBase(const TriaAccessorBase &)
void operator=(const TriaAccessorBase *)=delete
static const unsigned int structure_dimension
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
TriaAccessorBase(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *=nullptr)
void copy_from(const TriaAccessorBase &)
IteratorState::IteratorStates state() const
int index() const
bool operator<(const TriaAccessorBase &other) const
TriaAccessorBase & operator=(const TriaAccessorBase &)
::internal::TriangulationImplementation::TriaObjects & objects() const
int level() const
const Triangulation< dim, spacedim > * tria
const Triangulation< dim, spacedim > & get_triangulation() const
bool operator==(const TriaAccessorBase &) const
static unsigned int line_index(const unsigned int i)
static RefinementCase< 0 > refinement_case()
static unsigned int max_refinement_depth()
Point< spacedim > center() const
static unsigned int quad_index(const unsigned int i)
TriaAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
TriaAccessor(const Triangulation< 1, spacedim > *tria, const VertexKind vertex_kind, const unsigned int vertex_index)
static TriaIterator< TriaAccessor< 0, 1, spacedim > > isotropic_child(const unsigned int)
Return an invalid object.
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
unsigned int n_lines() const
bool operator!=(const TriaAccessor &) const
TriaAccessor(const Triangulation< 1, spacedim > *tria=nullptr, const int=0, const int=0, const AccessorData *=nullptr)
static unsigned int n_active_descendants()
static int isotropic_child_index(const unsigned int i)
Returns -1.
types::boundary_id boundary_id() const
static typename ::internal::TriangulationImplementation::Iterators< 1, spacedim >::line_iterator line(const unsigned int)
Point< spacedim > & vertex(const unsigned int i=0) const
static typename ::internal::TriangulationImplementation::Iterators< 1, spacedim >::quad_iterator quad(const unsigned int i)
ReferenceCell reference_cell() const
static unsigned int child_iterator_to_index(const TriaIterator< TriaAccessor< 0, 1, spacedim >> &)
Return an invalid unsigned integer.
unsigned int vertex_index(const unsigned int i=0) const
const Triangulation< 1, spacedim > & get_triangulation() const
void copy_from(const TriaAccessor &)
void set_all_boundary_ids(const types::boundary_id) const
static TriaIterator< TriaAccessor< 0, 1, spacedim > > child(const unsigned int)
Return an invalid object.
static IteratorState::IteratorStates state()
static int child_index(const unsigned int i)
Returns -1.
const Manifold< 1, spacedim > & get_manifold() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
types::manifold_id manifold_id() const
static unsigned int n_children()
TriaAccessor(const TriaAccessor< structdim2, dim2, spacedim2 > &)
void set_boundary_id(const types::boundary_id) const
unsigned int n_vertices() const
void set_manifold_id(const types::manifold_id)
void copy_from(const TriaAccessorBase< 0, 1, spacedim > &)
bool operator==(const TriaAccessor &) const
const Triangulation< 1, spacedim > * tria
double extent_in_direction(const unsigned int axis) const
TriaAccessor(const Triangulation< dim, spacedim > *tria, const unsigned int vertex_index)
static RefinementCase< 0 > refinement_case()
bool operator!=(const TriaAccessor &) const
static typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int)
void copy_from(const TriaAccessor &)
static typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i)
TriaAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
static int child_index(const unsigned int i)
Returns -1.
bool operator==(const TriaAccessor &) const
static unsigned int line_index(const unsigned int i)
unsigned int vertex_index(const unsigned int i=0) const
const Triangulation< dim, spacedim > * tria
static int isotropic_child_index(const unsigned int i)
Returns -1.
Point< spacedim > center(const bool respect_manifold=false, const bool interpolate_from_surrounding=false) const
TriaAccessor(const Triangulation< dim, spacedim > *tria=nullptr, const int level=0, const int index=0, const AccessorData *=nullptr)
static unsigned int n_children()
static TriaIterator< TriaAccessor< 0, dim, spacedim > > isotropic_child(const unsigned int)
Return an invalid object.
static unsigned int child_iterator_to_index(const TriaIterator< TriaAccessor< 0, dim, spacedim >> &)
Return an invalid unsigned integer.
static unsigned int max_refinement_depth()
static unsigned int quad_index(const unsigned int i)
const Triangulation< dim, spacedim > & get_triangulation() const
static unsigned int n_active_descendants()
static TriaIterator< TriaAccessor< 0, dim, spacedim > > child(const unsigned int)
Return an invalid object.
TriaAccessor(const TriaAccessor< structdim2, dim2, spacedim2 > &)
IteratorState::IteratorStates state() const
Point< spacedim > & vertex(const unsigned int i=0) const
void clear_children() const
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
void set_boundary_id_internal(const types::boundary_id id) const
void set_user_index(const unsigned int p) const
TriaAccessor(TriaAccessor &&)=default
void clear_user_pointer() const
unsigned int child_iterator_to_index(const TriaIterator< TriaAccessor< structdim, dim, spacedim >> &child) const
unsigned int n_active_descendants() const
void recursively_set_user_index(const unsigned int p) const
void clear_user_data() const
TriaAccessor(const TriaAccessor< structdim2, dim2, spacedim2 > &)
Point< spacedim > & vertex(const unsigned int i) const
Point< structdim > real_to_unit_cell_affine_approximation(const Point< spacedim > &point) const
void recursively_clear_user_index() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
unsigned int line_index(const unsigned int i) const
void recursively_set_user_pointer(void *p) const
unsigned int n_lines() const
double extent_in_direction(const unsigned int axis) const
Point< spacedim > intermediate_point(const Point< structdim > &coordinates) const
unsigned int n_vertices() const
bool has_children() const
void recursively_clear_user_flag() const
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int i) const
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
void set_bounding_object_indices(const std::initializer_list< unsigned int > &new_indices) const
void clear_user_flag() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
TriaAccessor & operator=(TriaAccessor &&)=default
TriaAccessor & operator=(const TriaAccessor &)=delete
unsigned int n_children() const
void recursively_set_user_flag() const
void set_boundary_id(const types::boundary_id) const
bool user_flag_set() const
void set_used_flag() const
types::manifold_id manifold_id() const
void set_user_flag() const
TriaAccessor(const TriaAccessor &)=default
TriaIterator< TriaAccessor< 0, dim, spacedim > > vertex_iterator(const unsigned int i) const
void clear_used_flag() const
std::pair< Point< spacedim >, double > enclosing_ball() const
Point< spacedim > barycenter() const
unsigned int vertex_index(const unsigned int i) const
int isotropic_child_index(const unsigned int i) const
Point< spacedim > center(const bool respect_manifold=false, const bool interpolate_from_surrounding=false) const
~TriaAccessor()=default
void set_refinement_case(const RefinementCase< structdim > &ref_case) const
void clear_user_index() const
double minimum_vertex_distance() const
double measure() const
void set_all_boundary_ids(const types::boundary_id) const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > isotropic_child(const unsigned int i) const
void set_bounding_object_indices(const std::initializer_list< int > &new_indices) const
unsigned int max_refinement_depth() const
void set_line_orientation(const unsigned int line, const bool orientation) const
bool is_translation_of(const TriaIterator< TriaAccessor< structdim, dim, spacedim >> &o) const
unsigned int quad_index(const unsigned int i) const
unsigned int user_index() const
int child_index(const unsigned int i) const
void set_user_pointer(void *p) const
BoundingBox< spacedim > bounding_box() const
void recursively_clear_user_pointer() const
ReferenceCell reference_cell() const
void * user_pointer() const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > child(const unsigned int i) const
TriaAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
const Manifold< dim, spacedim > & get_manifold() const
RefinementCase< structdim > refinement_case() const
TriaAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
bool used() const
void set_children(const unsigned int i, const int index) const
void clear_refinement_case() const
double diameter() const
bool at_boundary() const
types::boundary_id boundary_id() const
unsigned int n_faces() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_CXX20_REQUIRES(condition)
Definition: config.h:166
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
unsigned int level
Definition: grid_out.cc:4617
#define DeclException0(Exception0)
Definition: exceptions.h:467
static ::ExceptionBase & ExcRefineCellNotActive()
static ::ExceptionBase & ExcCantSetChildren(int arg1)
static ::ExceptionBase & ExcCellFlaggedForCoarsening()
static ::ExceptionBase & ExcCellNotActive()
static ::ExceptionBase & ExcCellHasNoParent()
static ::ExceptionBase & ExcNeighborIsNotCoarser()
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()
static ::ExceptionBase & ExcCellNotUsed()
#define Assert(cond, exc)
Definition: exceptions.h:1616
static ::ExceptionBase & ExcCellFlaggedForRefinement()
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:490
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:512
static ::ExceptionBase & ExcSetOnlyEvenChildren(int arg1)
static ::ExceptionBase & ExcCellHasNoChildren()
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcFacesHaveNoLevel()
static ::ExceptionBase & ExcNoPeriodicNeighbor()
static ::ExceptionBase & ExcNeighborIsCoarser()
void set_all_manifold_ids(const types::manifold_id) const
void set_all_manifold_ids(const types::manifold_id)
void set_manifold_id(const types::manifold_id) const
bool line_orientation(const unsigned int line) const
static bool face_orientation(const unsigned int face)
Always return false.
unsigned char combined_face_orientation(const unsigned int face) const
bool face_rotation(const unsigned int face) const
bool face_orientation(const unsigned int face) const
static bool face_flip(const unsigned int face)
Always return false.
static bool line_orientation(const unsigned int line)
Always return false.
static unsigned char combined_face_orientation(const unsigned int face)
Always return 0.
static bool face_flip(const unsigned int face)
Always return false.
static bool line_orientation(const unsigned int line)
Always return false.
static bool face_rotation(const unsigned int face)
Always return false.
static bool face_rotation(const unsigned int face)
Always return false.
static bool face_orientation(const unsigned int face)
Always return false.
static unsigned char combined_face_orientation(const unsigned int face)
Always return 0.
bool face_flip(const unsigned int face) const
void set_combined_face_orientation(const unsigned int face, const unsigned char combined_orientation) const
@ past_the_end
Iterator reached end of container.
@ valid
Iterator points to a valid object.
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:192
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
unsigned int manifold_id
Definition: types.h:153
unsigned int subdomain_id
Definition: types.h:44
unsigned int global_cell_index
Definition: types.h:118
unsigned int material_id
Definition: types.h:164
unsigned int boundary_id
Definition: types.h:141