Reference documentation for deal.II version GIT 1e9e64a91d 2022-09-28 19:20: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 - 2022 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_tria_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 
33 #include <boost/container/small_vector.hpp>
35 
36 #include <utility>
37 
38 
40 
41 // Forward declarations
42 #ifndef DOXYGEN
43 template <int dim, int spacedim>
44 class Triangulation;
45 template <typename Accessor>
46 class TriaRawIterator;
47 template <typename Accessor>
48 class TriaIterator;
49 template <typename Accessor>
50 class TriaActiveIterator;
51 
52 namespace parallel
53 {
54  template <int dim, int spacedim>
55  class TriangulationBase;
56 }
57 
58 template <int dim, int spacedim>
59 class DoFHandler;
60 template <int dim, int spacedim, bool lda>
61 class DoFCellAccessor;
62 
63 
64 template <int dim, int spacedim>
65 class Manifold;
66 
67 template <int dim, int spacedim>
68 class Mapping;
69 #endif
70 
71 namespace internal
72 {
73  namespace TriangulationImplementation
74  {
75  class TriaObjects;
76  struct Implementation;
77  struct ImplementationMixedMesh;
78  } // namespace TriangulationImplementation
79 
80  namespace TriaAccessorImplementation
81  {
82  struct Implementation;
83 
89  template <int structdim, int dim>
91  {
92  struct type
93  {
97  type() = default;
98 
102  type(const int level)
103  {
104  Assert(level == 0, ExcInternalError());
105  (void)level; // removes -Wunused-parameter warning in optimized mode
106  }
107 
111  operator int() const
112  {
113  return 0;
114  }
115 
116  void
117  operator++() const
118  {
119  Assert(false, ExcInternalError());
120  }
121 
122  void
123  operator--() const
124  {
125  Assert(false, ExcInternalError());
126  }
127  };
128  };
129 
130 
136  template <int dim>
137  struct PresentLevelType<dim, dim>
138  {
139  using type = int;
140  };
141  } // namespace TriaAccessorImplementation
142 } // namespace internal
143 template <int structdim, int dim, int spacedim>
144 class TriaAccessor;
145 template <int dim, int spacedim>
146 class TriaAccessor<0, dim, spacedim>;
147 template <int spacedim>
148 class TriaAccessor<0, 1, spacedim>;
149 
150 // note: the file tria_accessor.templates.h is included at the end of
151 // this file. this includes a lot of templates. originally, this was
152 // only done in debug mode, but led to cyclic reduction problems and
153 // so is now on by default.
154 
155 
160 {
165  "The operation you are attempting can only be performed for "
166  "(cell, face, or edge) iterators that point to valid "
167  "objects. These objects need not necessarily be active, "
168  "i.e., have no children, but they need to be part of a "
169  "triangulation. (The objects pointed to by an iterator "
170  "may -- after coarsening -- also be objects that used "
171  "to be part of a triangulation, but are now no longer "
172  "used. Their memory location may have been retained "
173  "for re-use upon the next mesh refinement, but is "
174  "currently unused.)");
185  "The operation you are attempting can only be performed for "
186  "(cell, face, or edge) iterators that point to 'active' "
187  "objects. 'Active' objects are those that do not have "
188  "children (in the case of cells), or that are part of "
189  "an active cell (in the case of faces or edges). However, "
190  "the object on which you are trying the current "
191  "operation is not 'active' in this sense.");
198  "The operation you are attempting can only be performed for "
199  "(cell, face, or edge) iterators that have children, "
200  "but the object on which you are trying the current "
201  "operation does not have any.");
209  "The operation you are attempting can only be performed for "
210  "(cell, face, or edge) iterators that have a parent object, "
211  "but the object on which you are trying the current "
212  "operation does not have one -- i.e., it is on the "
213  "coarsest level of the triangulation.");
218  int,
219  << "You can only set the child index if the cell does not "
220  << "currently have children registered; or you can clear it. "
221  << "The given index was " << arg1
222  << " (-1 means: clear children).");
226  template <typename AccessorType>
228  AccessorType,
229  << "You tried to dereference an iterator for which this "
230  << "is not possible. More information on this iterator: "
231  << "index=" << arg1.index() << ", state="
232  << (arg1.state() == IteratorState::valid ?
233  "valid" :
234  (arg1.state() == IteratorState::past_the_end ?
235  "past_the_end" :
236  "invalid")));
241  "Iterators can only be compared if they point to the same "
242  "triangulation, or if neither of them are associated "
243  "with a triangulation.");
244  // TODO: Write documentation!
249  // TODO: Write documentation!
269  // TODO: Write documentation!
275  int,
276  << "You can only set the child index of an even numbered child."
277  << "The number of the child given was " << arg1 << '.');
278 } // namespace TriaAccessorExceptions
279 
280 
306 template <int structdim, int dim, int spacedim = dim>
308 {
309 public:
315  static constexpr unsigned int space_dimension = spacedim;
316 
322  static constexpr unsigned int dimension = dim;
323 
329  static const unsigned int structure_dimension = structdim;
330 
340  void
341  operator=(const TriaAccessorBase *) = delete;
342 
343 protected:
349  using AccessorData = void;
350 
355  const int level = -1,
356  const int index = -1,
357  const AccessorData * = nullptr);
358 
363 
371  void
373 
379 
390  bool
391  operator<(const TriaAccessorBase &other) const;
392 
393 protected:
397  bool
398  operator==(const TriaAccessorBase &) const;
399 
403  bool
404  operator!=(const TriaAccessorBase &) const;
405 
419  void
421 
429  void
439  objects() const;
440 
441 public:
447  using LocalData = void *;
448 
472  int
473  level() const;
474 
501  int
502  index() const;
503 
509  state() const;
510 
517 
521 protected:
526  typename ::internal::TriaAccessorImplementation::
527  PresentLevelType<structdim, dim>::type present_level;
528 
534 
539 
540 private:
541  template <typename Accessor>
542  friend class TriaRawIterator;
543  template <typename Accessor>
544  friend class TriaIterator;
545  template <typename Accessor>
546  friend class TriaActiveIterator;
547 };
548 
549 
550 
571 template <int structdim, int dim, int spacedim = dim>
572 class InvalidAccessor : public TriaAccessorBase<structdim, dim, spacedim>
573 {
574 public:
578  using AccessorData =
580 
589  const int level = -1,
590  const int index = -1,
591  const AccessorData * local_data = nullptr);
592 
601 
606  template <typename OtherAccessor>
607  InvalidAccessor(const OtherAccessor &);
608 
612  void
614 
618  bool
619  operator==(const InvalidAccessor &) const;
620  bool
621  operator!=(const InvalidAccessor &) const;
622 
626  void
627  operator++() const;
628  void
629  operator--() const;
630 
635  bool
636  used() const;
637 
642  bool
643  has_children() const;
644 
649  manifold_id() const;
650 
654  unsigned int
655  user_index() const;
656 
660  void
661  set_user_index(const unsigned int p) const;
662 
666  void
668 
673  vertex(const unsigned int i) const;
674 
679  typename ::internal::TriangulationImplementation::
680  Iterators<dim, spacedim>::line_iterator
681  line(const unsigned int i) const;
682 
687  typename ::internal::TriangulationImplementation::
688  Iterators<dim, spacedim>::quad_iterator
689  quad(const unsigned int i) const;
690 };
691 
692 
693 
711 template <int structdim, int dim, int spacedim>
712 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
713 {
714 public:
718  using AccessorData =
720 
725  const int level = -1,
726  const int index = -1,
727  const AccessorData * local_data = nullptr);
728 
733  TriaAccessor(const TriaAccessor &) = default;
734 
738  TriaAccessor(TriaAccessor &&) = default; // NOLINT
739 
752  template <int structdim2, int dim2, int spacedim2>
754 
759  template <int structdim2, int dim2, int spacedim2>
761 
771  TriaAccessor &
772  operator=(const TriaAccessor &) = delete;
773 
777  TriaAccessor &
778  operator=(TriaAccessor &&) = default; // NOLINT
779 
783  ~TriaAccessor() = default;
784 
791  bool
792  used() const;
793 
806  vertex_iterator(const unsigned int i) const;
807 
823  unsigned int
824  vertex_index(const unsigned int i) const;
825 
864  vertex(const unsigned int i) const;
865 
869  typename ::internal::TriangulationImplementation::
870  Iterators<dim, spacedim>::line_iterator
871  line(const unsigned int i) const;
872 
879  unsigned int
880  line_index(const unsigned int i) const;
881 
885  typename ::internal::TriangulationImplementation::
886  Iterators<dim, spacedim>::quad_iterator
887  quad(const unsigned int i) const;
888 
895  unsigned int
896  quad_index(const unsigned int i) const;
912  unsigned char
913  combined_face_orientation(const unsigned int face) const;
914 
926  bool
927  face_orientation(const unsigned int face) const;
928 
938  bool
939  face_flip(const unsigned int face) const;
940 
950  bool
951  face_rotation(const unsigned int face) const;
952 
965  bool
966  line_orientation(const unsigned int line) const;
981  bool
982  has_children() const;
983 
988  unsigned int
989  n_children() const;
990 
995  unsigned int
997 
1011  unsigned int
1013 
1027  unsigned int
1029 
1034  child(const unsigned int i) const;
1035 
1040  unsigned int
1043 
1053  isotropic_child(const unsigned int i) const;
1054 
1060 
1066  int
1067  child_index(const unsigned int i) const;
1068 
1074  int
1075  isotropic_child_index(const unsigned int i) const;
1098  boundary_id() const;
1099 
1129  void
1131 
1162  void
1164 
1172  bool
1173  at_boundary() const;
1174 
1184  const Manifold<dim, spacedim> &
1185  get_manifold() const;
1186 
1208  manifold_id() const;
1209 
1227  void
1229 
1243  void
1245 
1262  bool
1263  user_flag_set() const;
1264 
1270  void
1271  set_user_flag() const;
1272 
1278  void
1280 
1286  void
1288 
1294  void
1296 
1302  void
1304 
1316  void
1317  set_user_pointer(void *p) const;
1318 
1324  void
1326 
1342  void *
1343  user_pointer() const;
1344 
1366  void
1368 
1375  void
1377 
1387  void
1388  set_user_index(const unsigned int p) const;
1389 
1395  void
1397 
1409  unsigned int
1410  user_index() const;
1411 
1429  void
1430  recursively_set_user_index(const unsigned int p) const;
1431 
1440  void
1476  double
1477  diameter() const;
1478 
1505  std::pair<Point<spacedim>, double>
1507 
1517  bounding_box() const;
1518 
1528  double
1529  extent_in_direction(const unsigned int axis) const;
1530 
1534  double
1536 
1551  intermediate_point(const Point<structdim> &coordinates) const;
1552 
1577 
1613  center(const bool respect_manifold = false,
1614  const bool interpolate_from_surrounding = false) const;
1615 
1634  barycenter() const;
1635 
1661  double
1662  measure() const;
1663 
1678  bool
1681 
1687 
1691  unsigned int
1692  n_vertices() const;
1693 
1697  unsigned int
1698  n_lines() const;
1699 
1705  unsigned int
1706  n_faces() const;
1707 
1714 
1720  line_indices() const;
1721 
1729  face_indices() const;
1730 
1735 private:
1740  void
1742 
1750  void
1752  const std::initializer_list<int> &new_indices) const;
1753 
1757  void
1759  const std::initializer_list<unsigned int> &new_indices) const;
1760 
1768  void
1769  set_line_orientation(const unsigned int line, const bool orientation) const;
1770 
1781  void
1782  set_face_orientation(const unsigned int face, const bool orientation) const;
1783 
1790  void
1791  set_face_flip(const unsigned int face, const bool flip) const;
1792 
1799  void
1800  set_face_rotation(const unsigned int face, const bool rotation) const;
1801 
1805  void
1806  set_used_flag() const;
1807 
1811  void
1813 
1822  void
1824 
1832  void
1834 
1841  void
1842  set_children(const unsigned int i, const int index) const;
1843 
1848  void
1850 
1851 private:
1852  template <int, int>
1853  friend class Triangulation;
1854 
1855  friend struct ::internal::TriangulationImplementation::Implementation;
1856  friend struct ::internal::TriangulationImplementation::
1857  ImplementationMixedMesh;
1858  friend struct ::internal::TriaAccessorImplementation::Implementation;
1859 };
1860 
1861 
1862 
1881 template <int dim, int spacedim>
1882 class TriaAccessor<0, dim, spacedim>
1883 {
1884 public:
1890  static constexpr unsigned int space_dimension = spacedim;
1891 
1897  static constexpr unsigned int dimension = dim;
1898 
1904  static const unsigned int structure_dimension = 0;
1905 
1909  using AccessorData = void;
1910 
1916  const unsigned int vertex_index);
1917 
1924  const int level = 0,
1925  const int index = 0,
1926  const AccessorData * = nullptr);
1927 
1931  template <int structdim2, int dim2, int spacedim2>
1933 
1937  template <int structdim2, int dim2, int spacedim2>
1939 
1944  state() const;
1945 
1950  static int
1952 
1957  int
1958  index() const;
1959 
1966 
1976  void
1978 
1982  void
1987  bool
1988  operator==(const TriaAccessor &) const;
1989 
1993  bool
1994  operator!=(const TriaAccessor &) const;
1995 
2023  unsigned int
2024  vertex_index(const unsigned int i = 0) const;
2025 
2031  Point<spacedim> &
2032  vertex(const unsigned int i = 0) const;
2033 
2038  typename ::internal::TriangulationImplementation::
2039  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2040 
2044  static unsigned int
2045  line_index(const unsigned int i);
2046 
2050  static typename ::internal::TriangulationImplementation::
2051  Iterators<dim, spacedim>::quad_iterator
2052  quad(const unsigned int i);
2053 
2057  static unsigned int
2058  quad_index(const unsigned int i);
2059 
2075  double
2076  diameter() const;
2077 
2085  double
2086  extent_in_direction(const unsigned int axis) const;
2087 
2096  center(const bool respect_manifold = false,
2097  const bool interpolate_from_surrounding = false) const;
2098 
2107  double
2108  measure() const;
2123  static unsigned char
2124  combined_face_orientation(const unsigned int face);
2125 
2129  static bool
2130  face_orientation(const unsigned int face);
2131 
2135  static bool
2136  face_flip(const unsigned int face);
2137 
2141  static bool
2142  face_rotation(const unsigned int face);
2143 
2147  static bool
2148  line_orientation(const unsigned int line);
2149 
2164  static bool
2166 
2171  static unsigned int
2173 
2178  static unsigned int
2180 
2185  static unsigned int
2187 
2188 
2192  static unsigned int
2194 
2198  static unsigned int
2200 
2205  child(const unsigned int);
2206 
2211  isotropic_child(const unsigned int);
2212 
2216  static RefinementCase<0>
2218 
2222  static int
2223  child_index(const unsigned int i);
2224 
2228  static int
2229  isotropic_child_index(const unsigned int i);
2237  bool
2238  used() const;
2239 
2240 protected:
2248  void
2250 
2259  bool
2260  operator<(const TriaAccessor &other) const;
2261 
2266 
2270  unsigned int global_vertex_index;
2271 
2272 private:
2273  template <typename Accessor>
2274  friend class TriaRawIterator;
2275  template <typename Accessor>
2276  friend class TriaIterator;
2277  template <typename Accessor>
2278  friend class TriaActiveIterator;
2279 };
2280 
2281 
2282 
2299 template <int spacedim>
2300 class TriaAccessor<0, 1, spacedim>
2301 {
2302 public:
2308  static constexpr unsigned int space_dimension = spacedim;
2309 
2315  static constexpr unsigned int dimension = 1;
2316 
2322  static const unsigned int structure_dimension = 0;
2323 
2327  using AccessorData = void;
2328 
2334  {
2346  right_vertex
2347  };
2348 
2361  const VertexKind vertex_kind,
2362  const unsigned int vertex_index);
2363 
2370  const int = 0,
2371  const int = 0,
2372  const AccessorData * = nullptr);
2373 
2377  template <int structdim2, int dim2, int spacedim2>
2379 
2383  template <int structdim2, int dim2, int spacedim2>
2385 
2390  void
2392 
2398  void
2400 
2408 
2413  static int
2415 
2420  int
2421  index() const;
2422 
2429 
2440  void
2441  operator++() const;
2442 
2447  void
2448  operator--() const;
2452  bool
2453  operator==(const TriaAccessor &) const;
2454 
2458  bool
2459  operator!=(const TriaAccessor &) const;
2460 
2469  bool
2470  operator<(const TriaAccessor &other) const;
2471 
2498  unsigned int
2499  vertex_index(const unsigned int i = 0) const;
2500 
2506  Point<spacedim> &
2507  vertex(const unsigned int i = 0) const;
2508 
2514  center() const;
2515 
2520  typename ::internal::TriangulationImplementation::
2521  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2522 
2529  static unsigned int
2530  line_index(const unsigned int i);
2531 
2535  static typename ::internal::TriangulationImplementation::
2536  Iterators<1, spacedim>::quad_iterator
2537  quad(const unsigned int i);
2538 
2545  static unsigned int
2546  quad_index(const unsigned int i);
2547 
2557  bool
2558  at_boundary() const;
2559 
2575  boundary_id() const;
2576 
2580  const Manifold<1, spacedim> &
2581  get_manifold() const;
2582 
2590  manifold_id() const;
2591 
2592 
2604  bool
2605  user_flag_set() const;
2606 
2612  void
2613  set_user_flag() const;
2614 
2620  void
2621  clear_user_flag() const;
2622 
2628  void
2629  recursively_set_user_flag() const;
2630 
2636  void
2638 
2644  void
2645  clear_user_data() const;
2646 
2658  void
2659  set_user_pointer(void *p) const;
2660 
2666  void
2667  clear_user_pointer() const;
2668 
2684  void *
2685  user_pointer() const;
2686 
2708  void
2709  recursively_set_user_pointer(void *p) const;
2710 
2717  void
2719 
2729  void
2730  set_user_index(const unsigned int p) const;
2731 
2737  void
2738  clear_user_index() const;
2739 
2751  unsigned int
2752  user_index() const;
2753 
2771  void
2772  recursively_set_user_index(const unsigned int p) const;
2773 
2782  void
2798  static unsigned char
2799  combined_face_orientation(const unsigned int face);
2800 
2804  static bool
2805  face_orientation(const unsigned int face);
2806 
2810  static bool
2811  face_flip(const unsigned int face);
2812 
2816  static bool
2817  face_rotation(const unsigned int face);
2818 
2822  static bool
2823  line_orientation(const unsigned int line);
2824 
2839  static bool
2841 
2846  static unsigned int
2848 
2853  static unsigned int
2855 
2860  static unsigned int
2862 
2863 
2867  static unsigned int
2869 
2873  static unsigned int
2875 
2880  child(const unsigned int);
2881 
2886  isotropic_child(const unsigned int);
2887 
2891  static RefinementCase<0>
2893 
2897  static int
2898  child_index(const unsigned int i);
2899 
2903  static int
2904  isotropic_child_index(const unsigned int i);
2937  void
2939 
2946  void
2948 
2960  void
2962 
2974  void
2983  bool
2984  used() const;
2985 
2991 
2995  unsigned int
2996  n_vertices() const;
2997 
3001  unsigned int
3002  n_lines() const;
3003 
3010 
3016  line_indices() const;
3017 
3018 protected:
3023 
3029 
3033  unsigned int global_vertex_index;
3034 };
3035 
3036 
3037 
3053 template <int dim, int spacedim = dim>
3054 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
3055 {
3056 public:
3061 
3066 
3078  const int level = -1,
3079  const int index = -1,
3080  const AccessorData * local_data = nullptr);
3081 
3086 
3099  template <int structdim2, int dim2, int spacedim2>
3101 
3106  template <int structdim2, int dim2, int spacedim2>
3108 
3113 
3117  // NOLINTNEXTLINE OSX does not compile with noexcept
3119 
3123  ~CellAccessor() = default;
3124 
3136 
3140  // NOLINTNEXTLINE OSX does not compile with noexcept
3142  operator=(CellAccessor<dim, spacedim> &&) = default; // NOLINT
3143 
3167  as_dof_handler_iterator(const DoFHandler<dim, spacedim> &dof_handler) const;
3168 
3185  child(const unsigned int i) const;
3186 
3190  boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
3193 
3197  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
3198  face(const unsigned int i) const;
3199 
3204  unsigned int
3207 
3211  boost::container::small_vector<
3212  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
3215 
3225  unsigned int
3226  face_index(const unsigned int i) const;
3227 
3276  neighbor_child_on_subface(const unsigned int face_no,
3277  const unsigned int subface_no) const;
3278 
3304  neighbor(const unsigned int face_no) const;
3305 
3313  int
3314  neighbor_index(const unsigned int face_no) const;
3315 
3323  int
3324  neighbor_level(const unsigned int face_no) const;
3325 
3337  unsigned int
3338  neighbor_of_neighbor(const unsigned int face_no) const;
3339 
3350  bool
3351  neighbor_is_coarser(const unsigned int face_no) const;
3352 
3367  std::pair<unsigned int, unsigned int>
3368  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3369 
3376  unsigned int
3377  neighbor_face_no(const unsigned int neighbor) const;
3378 
3382  static bool
3384 
3398  bool
3399  has_periodic_neighbor(const unsigned int i) const;
3400 
3418  periodic_neighbor(const unsigned int i) const;
3419 
3428  neighbor_or_periodic_neighbor(const unsigned int i) const;
3429 
3445  periodic_neighbor_child_on_subface(const unsigned int face_no,
3446  const unsigned int subface_no) const;
3447 
3458  std::pair<unsigned int, unsigned int>
3459  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3460 
3466  int
3467  periodic_neighbor_index(const unsigned int i) const;
3468 
3474  int
3475  periodic_neighbor_level(const unsigned int i) const;
3476 
3491  unsigned int
3492  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3493 
3499  unsigned int
3500  periodic_neighbor_face_no(const unsigned int i) const;
3501 
3508  bool
3509  periodic_neighbor_is_coarser(const unsigned int i) const;
3510 
3527  bool
3528  at_boundary(const unsigned int i) const;
3529 
3538  bool
3539  at_boundary() const;
3540 
3548  bool
3549  has_boundary_lines() const;
3577 
3595  void
3598 
3602  void
3604 
3612  bool
3614  const unsigned int face_no,
3615  const RefinementCase<dim - 1> &face_refinement_case =
3617 
3623  bool
3624  flag_for_line_refinement(const unsigned int line_no) const;
3625 
3635  subface_case(const unsigned int face_no) const;
3636 
3640  bool
3642 
3647  void
3649 
3653  void
3678  material_id() const;
3679 
3691  void
3692  set_material_id(const types::material_id new_material_id) const;
3693 
3702  void
3703  recursively_set_material_id(const types::material_id new_material_id) const;
3730  subdomain_id() const;
3731 
3747  void
3748  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3749 
3760 
3765  void
3767  const types::subdomain_id new_level_subdomain_id) const;
3768 
3769 
3785  void
3787  const types::subdomain_id new_subdomain_id) const;
3811 
3821 
3835  bool
3836  direction_flag() const;
3837 
3863  unsigned int
3865 
3873  int
3874  parent_index() const;
3875 
3882  parent() const;
3883 
3903  bool
3904  is_active() const;
3905 
3925  bool
3927 
3932  bool
3934 
3968  bool
3969  is_ghost() const;
3970 
3976  bool
3978 
4005  bool
4006  is_artificial() const;
4007 
4014  bool
4016 
4030  bool
4032 
4041  void
4042  set_neighbor(const unsigned int i,
4043  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
4044 
4058  CellId
4059  id() const;
4060 
4062 
4066  double
4067  diameter(const Mapping<dim, spacedim> &mapping) const;
4068 
4086 
4087 protected:
4103  unsigned int
4104  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
4105 
4111  template <int dim_, int spacedim_>
4112  bool
4113  point_inside_codim(const Point<spacedim_> &p) const;
4114 
4115 
4116 
4117 private:
4122  void
4123  set_active_cell_index(const unsigned int active_cell_index) const;
4124 
4128  void
4130 
4134  void
4136 
4140  void
4141  set_parent(const unsigned int parent_index);
4142 
4149  void
4150  set_direction_flag(const bool new_direction_flag) const;
4151 
4152  template <int, int>
4153  friend class Triangulation;
4154 
4155  template <int, int>
4157 
4158  friend struct ::internal::TriangulationImplementation::Implementation;
4159  friend struct ::internal::TriangulationImplementation::
4160  ImplementationMixedMesh;
4161 };
4162 
4163 
4164 
4165 /* ----- declaration of explicit specializations and general templates ----- */
4166 
4167 
4168 template <int structdim, int dim, int spacedim>
4169 template <typename OtherAccessor>
4171  const OtherAccessor &)
4172 {
4173  Assert(false,
4174  ExcMessage("You are attempting an illegal conversion between "
4175  "iterator/accessor types. The constructor you call "
4176  "only exists to make certain template constructs "
4177  "easier to write as dimension independent code but "
4178  "the conversion is not valid in the current context."));
4179 }
4180 
4181 
4182 
4183 template <int structdim, int dim, int spacedim>
4184 template <int structdim2, int dim2, int spacedim2>
4187 {
4188  Assert(false,
4189  ExcMessage("You are attempting an illegal conversion between "
4190  "iterator/accessor types. The constructor you call "
4191  "only exists to make certain template constructs "
4192  "easier to write as dimension independent code but "
4193  "the conversion is not valid in the current context."));
4194 }
4195 
4196 
4197 
4198 template <int dim, int spacedim>
4199 template <int structdim2, int dim2, int spacedim2>
4202 {
4203  Assert(false,
4204  ExcMessage("You are attempting an illegal conversion between "
4205  "iterator/accessor types. The constructor you call "
4206  "only exists to make certain template constructs "
4207  "easier to write as dimension independent code but "
4208  "the conversion is not valid in the current context."));
4209 }
4210 
4211 
4212 
4213 template <int structdim, int dim, int spacedim>
4214 template <int structdim2, int dim2, int spacedim2>
4217 {
4218  Assert(false,
4219  ExcMessage("You are attempting an illegal conversion between "
4220  "iterator/accessor types. The constructor you call "
4221  "only exists to make certain template constructs "
4222  "easier to write as dimension independent code but "
4223  "the conversion is not valid in the current context."));
4224 }
4225 
4226 
4227 
4228 template <int dim, int spacedim>
4229 template <int structdim2, int dim2, int spacedim2>
4232 {
4233  Assert(false,
4234  ExcMessage("You are attempting an illegal conversion between "
4235  "iterator/accessor types. The constructor you call "
4236  "only exists to make certain template constructs "
4237  "easier to write as dimension independent code but "
4238  "the conversion is not valid in the current context."));
4239 }
4240 
4241 
4242 #ifndef DOXYGEN
4243 
4244 template <>
4245 bool
4247 template <>
4248 bool
4250 template <>
4251 bool
4253 template <>
4254 bool
4256 template <>
4257 bool
4259 template <>
4260 bool
4262 // -------------------------------------------------------------------
4263 
4264 template <>
4265 void
4267 
4268 #endif // DOXYGEN
4269 
4271 
4272 // include more templates in debug and optimized mode
4273 #include "tria_accessor.templates.h"
4274 
4275 #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
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:71
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
void set_manifold_id(const types::manifold_id) const
InvalidAccessor(const InvalidAccessor &)
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int i) const
InvalidAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
bool operator==(const InvalidAccessor &) const
void operator++() const
bool operator!=(const InvalidAccessor &) const
bool has_children() const
unsigned int user_index() const
Point< spacedim > & vertex(const unsigned int i) const
void copy_from(const InvalidAccessor &)
void set_user_index(const unsigned int p) const
types::manifold_id manifold_id() const
void operator--() const
bool used() const
InvalidAccessor(const OtherAccessor &)
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
Abstract base class for mapping classes.
Definition: mapping.h:311
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 bool face_orientation(const unsigned int face)
Always return false.
static unsigned int number_of_children()
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 > &)
static bool face_flip(const unsigned int face)
Always return false.
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 bool line_orientation(const unsigned int line)
Always return false.
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 bool face_rotation(const unsigned int face)
Always return false.
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 > &)
static unsigned char combined_face_orientation(const unsigned int face)
Always return 0.
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 char combined_face_orientation(const unsigned int face)
Always return 0.
static unsigned int n_children()
static bool face_flip(const unsigned int face)
Always return false.
static TriaIterator< TriaAccessor< 0, dim, spacedim > > isotropic_child(const unsigned int)
Return an invalid object.
static bool line_orientation(const unsigned int line)
Always return false.
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)
static bool face_rotation(const unsigned int face)
Always return false.
const Triangulation< dim, spacedim > & get_triangulation() const
static bool face_orientation(const unsigned int face)
Always return false.
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
static unsigned int number_of_children()
void clear_children() const
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
bool line_orientation(const unsigned int line) const
void set_boundary_id_internal(const types::boundary_id id) const
void set_user_index(const unsigned int p) const
void set_face_flip(const unsigned int face, const bool flip) 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
unsigned char combined_face_orientation(const unsigned int face) const
void set_face_orientation(const unsigned int face, const bool orientation) const
void set_face_rotation(const unsigned int face, const bool rotation) const
bool face_rotation(const unsigned int face) 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
bool face_orientation(const unsigned int face) 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
Point< spacedim > barycenter() const
BoundingBox< spacedim > bounding_box() 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
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 number_of_children() 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
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
bool face_flip(const unsigned int face) 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_DEPRECATED
Definition: config.h:164
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:457
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
Definition: config.h:471
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:458
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
Definition: config.h:510
unsigned int level
Definition: grid_out.cc:4608
#define DeclException0(Exception0)
Definition: exceptions.h:464
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:1501
static ::ExceptionBase & ExcCellFlaggedForRefinement()
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
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
@ 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:189
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
unsigned int manifold_id
Definition: types.h:146
unsigned int subdomain_id
Definition: types.h:43
unsigned int global_cell_index
Definition: types.h:110
unsigned int material_id
Definition: types.h:157
unsigned int boundary_id
Definition: types.h:134