Reference documentation for deal.II version GIT 7deb6c54a6 2023-06-09 18: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 - 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 
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  class TriangulationBase;
55 }
56 
57 template <int dim, int spacedim>
58 DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, 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>
573 {
574 public:
580  static constexpr unsigned int space_dimension = spacedim;
581 
587  static constexpr unsigned int dimension = dim;
588 
594  static const unsigned int structure_dimension = structdim;
595 
601  using AccessorData = void;
602 
610  InvalidAccessor(const void * parent = nullptr,
611  const int level = -1,
612  const int index = -1,
613  const AccessorData *local_data = nullptr);
614 
623 
628  template <typename OtherAccessor>
629  InvalidAccessor(const OtherAccessor &);
630 
634  void
636 
640  bool
641  operator==(const InvalidAccessor &) const;
642  bool
643  operator!=(const InvalidAccessor &) const;
644 
648  void
649  operator++() const;
650  void
651  operator--() const;
652 
658  state();
659 
660 
665  static int
666  level();
667 
672  static int
673  index();
674 
679  bool
680  used() const;
681 
686  bool
687  has_children() const;
688 
693  manifold_id() const;
694 
698  unsigned int
699  user_index() const;
700 
704  void
705  set_user_index(const unsigned int p) const;
706 
710  void
712 
717  vertex(const unsigned int i) const;
718 
723  void *
724  line(const unsigned int i) const;
725 
730  void *
731  quad(const unsigned int i) const;
732 };
733 
734 
735 
753 template <int structdim, int dim, int spacedim>
754 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
755 {
756 public:
760  using AccessorData =
762 
767  const int level = -1,
768  const int index = -1,
769  const AccessorData * local_data = nullptr);
770 
775  TriaAccessor(const TriaAccessor &) = default;
776 
780  TriaAccessor(TriaAccessor &&) = default; // NOLINT
781 
794  template <int structdim2, int dim2, int spacedim2>
796 
801  template <int structdim2, int dim2, int spacedim2>
803 
813  TriaAccessor &
814  operator=(const TriaAccessor &) = delete;
815 
819  TriaAccessor &
820  operator=(TriaAccessor &&) = default; // NOLINT
821 
825  ~TriaAccessor() = default;
826 
833  bool
834  used() const;
835 
848  vertex_iterator(const unsigned int i) const;
849 
865  unsigned int
866  vertex_index(const unsigned int i) const;
867 
906  vertex(const unsigned int i) const;
907 
911  typename ::internal::TriangulationImplementation::
912  Iterators<dim, spacedim>::line_iterator
913  line(const unsigned int i) const;
914 
921  unsigned int
922  line_index(const unsigned int i) const;
923 
927  typename ::internal::TriangulationImplementation::
928  Iterators<dim, spacedim>::quad_iterator
929  quad(const unsigned int i) const;
930 
937  unsigned int
938  quad_index(const unsigned int i) const;
954  unsigned char
955  combined_face_orientation(const unsigned int face) const;
956 
968  bool
969  face_orientation(const unsigned int face) const;
970 
980  bool
981  face_flip(const unsigned int face) const;
982 
992  bool
993  face_rotation(const unsigned int face) const;
994 
1007  bool
1008  line_orientation(const unsigned int line) const;
1023  bool
1024  has_children() const;
1025 
1030  unsigned int
1031  n_children() const;
1032 
1037  unsigned int
1039 
1053  unsigned int
1055 
1069  unsigned int
1071 
1076  child(const unsigned int i) const;
1077 
1082  unsigned int
1085 
1095  isotropic_child(const unsigned int i) const;
1096 
1102 
1108  int
1109  child_index(const unsigned int i) const;
1110 
1116  int
1117  isotropic_child_index(const unsigned int i) const;
1140  boundary_id() const;
1141 
1171  void
1173 
1204  void
1206 
1214  bool
1215  at_boundary() const;
1216 
1226  const Manifold<dim, spacedim> &
1227  get_manifold() const;
1228 
1250  manifold_id() const;
1251 
1269  void
1271 
1285  void
1287 
1304  bool
1305  user_flag_set() const;
1306 
1312  void
1313  set_user_flag() const;
1314 
1320  void
1322 
1328  void
1330 
1336  void
1338 
1344  void
1346 
1358  void
1359  set_user_pointer(void *p) const;
1360 
1366  void
1368 
1384  void *
1385  user_pointer() const;
1386 
1408  void
1410 
1417  void
1419 
1429  void
1430  set_user_index(const unsigned int p) const;
1431 
1437  void
1439 
1451  unsigned int
1452  user_index() const;
1453 
1471  void
1472  recursively_set_user_index(const unsigned int p) const;
1473 
1482  void
1518  double
1519  diameter() const;
1520 
1547  std::pair<Point<spacedim>, double>
1549 
1559  bounding_box() const;
1560 
1570  double
1571  extent_in_direction(const unsigned int axis) const;
1572 
1576  double
1578 
1593  intermediate_point(const Point<structdim> &coordinates) const;
1594 
1619 
1655  center(const bool respect_manifold = false,
1656  const bool interpolate_from_surrounding = false) const;
1657 
1676  barycenter() const;
1677 
1703  double
1704  measure() const;
1705 
1720  bool
1723 
1729 
1733  unsigned int
1734  n_vertices() const;
1735 
1739  unsigned int
1740  n_lines() const;
1741 
1747  unsigned int
1748  n_faces() const;
1749 
1756 
1762  line_indices() const;
1763 
1771  face_indices() const;
1772 
1777 private:
1782  void
1784 
1792  void
1794  const std::initializer_list<int> &new_indices) const;
1795 
1799  void
1801  const std::initializer_list<unsigned int> &new_indices) const;
1802 
1810  void
1811  set_line_orientation(const unsigned int line, const bool orientation) const;
1812 
1820  void
1821  set_combined_face_orientation(const unsigned int face,
1822  const unsigned char combined_orientation) const;
1823 
1827  void
1828  set_used_flag() const;
1829 
1833  void
1835 
1844  void
1846 
1854  void
1856 
1863  void
1864  set_children(const unsigned int i, const int index) const;
1865 
1870  void
1872 
1873 private:
1874  friend class Triangulation<dim, spacedim>;
1875 
1876  friend struct ::internal::TriangulationImplementation::Implementation;
1877  friend struct ::internal::TriangulationImplementation::
1878  ImplementationMixedMesh;
1879  friend struct ::internal::TriaAccessorImplementation::Implementation;
1880 };
1881 
1882 
1883 
1902 template <int dim, int spacedim>
1903 class TriaAccessor<0, dim, spacedim>
1904 {
1905 public:
1911  static constexpr unsigned int space_dimension = spacedim;
1912 
1918  static constexpr unsigned int dimension = dim;
1919 
1925  static const unsigned int structure_dimension = 0;
1926 
1930  using AccessorData = void;
1931 
1937  const unsigned int vertex_index);
1938 
1945  const int level = 0,
1946  const int index = 0,
1947  const AccessorData * = nullptr);
1948 
1952  template <int structdim2, int dim2, int spacedim2>
1954 
1958  template <int structdim2, int dim2, int spacedim2>
1960 
1965  state() const;
1966 
1971  static int
1973 
1978  int
1979  index() const;
1980 
1987 
1997  void
1999 
2003  void
2008  bool
2009  operator==(const TriaAccessor &) const;
2010 
2014  bool
2015  operator!=(const TriaAccessor &) const;
2016 
2044  unsigned int
2045  vertex_index(const unsigned int i = 0) const;
2046 
2052  Point<spacedim> &
2053  vertex(const unsigned int i = 0) const;
2054 
2059  typename ::internal::TriangulationImplementation::
2060  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2061 
2065  static unsigned int
2066  line_index(const unsigned int i);
2067 
2071  static typename ::internal::TriangulationImplementation::
2072  Iterators<dim, spacedim>::quad_iterator
2073  quad(const unsigned int i);
2074 
2078  static unsigned int
2079  quad_index(const unsigned int i);
2080 
2096  double
2097  diameter() const;
2098 
2106  double
2107  extent_in_direction(const unsigned int axis) const;
2108 
2117  center(const bool respect_manifold = false,
2118  const bool interpolate_from_surrounding = false) const;
2119 
2128  double
2129  measure() const;
2144  static unsigned char
2145  combined_face_orientation(const unsigned int face);
2146 
2150  static bool
2151  face_orientation(const unsigned int face);
2152 
2156  static bool
2157  face_flip(const unsigned int face);
2158 
2162  static bool
2163  face_rotation(const unsigned int face);
2164 
2168  static bool
2169  line_orientation(const unsigned int line);
2170 
2185  static bool
2187 
2192  static unsigned int
2194 
2199  static unsigned int
2201 
2206  static unsigned int
2208 
2209 
2213  static unsigned int
2215 
2219  static unsigned int
2221 
2226  child(const unsigned int);
2227 
2232  isotropic_child(const unsigned int);
2233 
2237  static RefinementCase<0>
2239 
2243  static int
2244  child_index(const unsigned int i);
2245 
2249  static int
2250  isotropic_child_index(const unsigned int i);
2258  bool
2259  used() const;
2260 
2261 protected:
2269  void
2271 
2280  bool
2281  operator<(const TriaAccessor &other) const;
2282 
2287 
2291  unsigned int global_vertex_index;
2292 
2293 private:
2294  template <typename Accessor>
2295  friend class TriaRawIterator;
2296  template <typename Accessor>
2297  friend class TriaIterator;
2298  template <typename Accessor>
2299  friend class TriaActiveIterator;
2300 };
2301 
2302 
2303 
2320 template <int spacedim>
2321 class TriaAccessor<0, 1, spacedim>
2322 {
2323 public:
2329  static constexpr unsigned int space_dimension = spacedim;
2330 
2336  static constexpr unsigned int dimension = 1;
2337 
2343  static const unsigned int structure_dimension = 0;
2344 
2348  using AccessorData = void;
2349 
2355  {
2367  right_vertex
2368  };
2369 
2382  const VertexKind vertex_kind,
2383  const unsigned int vertex_index);
2384 
2391  const int = 0,
2392  const int = 0,
2393  const AccessorData * = nullptr);
2394 
2398  template <int structdim2, int dim2, int spacedim2>
2400 
2404  template <int structdim2, int dim2, int spacedim2>
2406 
2411  void
2413 
2419  void
2421 
2429 
2434  static int
2436 
2441  int
2442  index() const;
2443 
2450 
2461  void
2462  operator++() const;
2463 
2468  void
2469  operator--() const;
2473  bool
2474  operator==(const TriaAccessor &) const;
2475 
2479  bool
2480  operator!=(const TriaAccessor &) const;
2481 
2490  bool
2491  operator<(const TriaAccessor &other) const;
2492 
2519  unsigned int
2520  vertex_index(const unsigned int i = 0) const;
2521 
2527  Point<spacedim> &
2528  vertex(const unsigned int i = 0) const;
2529 
2535  center() const;
2536 
2541  typename ::internal::TriangulationImplementation::
2542  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2543 
2550  static unsigned int
2551  line_index(const unsigned int i);
2552 
2556  static typename ::internal::TriangulationImplementation::
2557  Iterators<1, spacedim>::quad_iterator
2558  quad(const unsigned int i);
2559 
2566  static unsigned int
2567  quad_index(const unsigned int i);
2568 
2578  bool
2579  at_boundary() const;
2580 
2596  boundary_id() const;
2597 
2601  const Manifold<1, spacedim> &
2602  get_manifold() const;
2603 
2611  manifold_id() const;
2612 
2613 
2625  bool
2626  user_flag_set() const;
2627 
2633  void
2634  set_user_flag() const;
2635 
2641  void
2642  clear_user_flag() const;
2643 
2649  void
2650  recursively_set_user_flag() const;
2651 
2657  void
2659 
2665  void
2666  clear_user_data() const;
2667 
2679  void
2680  set_user_pointer(void *p) const;
2681 
2687  void
2688  clear_user_pointer() const;
2689 
2705  void *
2706  user_pointer() const;
2707 
2729  void
2730  recursively_set_user_pointer(void *p) const;
2731 
2738  void
2740 
2750  void
2751  set_user_index(const unsigned int p) const;
2752 
2758  void
2759  clear_user_index() const;
2760 
2772  unsigned int
2773  user_index() const;
2774 
2792  void
2793  recursively_set_user_index(const unsigned int p) const;
2794 
2803  void
2819  static unsigned char
2820  combined_face_orientation(const unsigned int face);
2821 
2825  static bool
2826  face_orientation(const unsigned int face);
2827 
2831  static bool
2832  face_flip(const unsigned int face);
2833 
2837  static bool
2838  face_rotation(const unsigned int face);
2839 
2843  static bool
2844  line_orientation(const unsigned int line);
2845 
2860  static bool
2862 
2867  static unsigned int
2869 
2874  static unsigned int
2876 
2881  static unsigned int
2883 
2884 
2888  static unsigned int
2890 
2894  static unsigned int
2896 
2901  child(const unsigned int);
2902 
2907  isotropic_child(const unsigned int);
2908 
2912  static RefinementCase<0>
2914 
2918  static int
2919  child_index(const unsigned int i);
2920 
2924  static int
2925  isotropic_child_index(const unsigned int i);
2958  void
2960 
2967  void
2969 
2981  void
2983 
2995  void
3004  bool
3005  used() const;
3006 
3012 
3016  unsigned int
3017  n_vertices() const;
3018 
3022  unsigned int
3023  n_lines() const;
3024 
3031 
3037  line_indices() const;
3038 
3039 protected:
3044 
3050 
3054  unsigned int global_vertex_index;
3055 };
3056 
3057 
3058 
3074 template <int dim, int spacedim = dim>
3075 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
3076 {
3077 public:
3082 
3087 
3099  const int level = -1,
3100  const int index = -1,
3101  const AccessorData * local_data = nullptr);
3102 
3107 
3120  template <int structdim2, int dim2, int spacedim2>
3122 
3127  template <int structdim2, int dim2, int spacedim2>
3129 
3134 
3138  // NOLINTNEXTLINE OSX does not compile with noexcept
3140 
3144  ~CellAccessor() = default;
3145 
3157 
3161  // NOLINTNEXTLINE OSX does not compile with noexcept
3163  operator=(CellAccessor<dim, spacedim> &&) = default; // NOLINT
3164 
3188  as_dof_handler_iterator(const DoFHandler<dim, spacedim> &dof_handler) const;
3189 
3206  child(const unsigned int i) const;
3207 
3211  boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
3214 
3218  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
3219  face(const unsigned int i) const;
3220 
3225  unsigned int
3228 
3232  boost::container::small_vector<
3233  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
3236 
3246  unsigned int
3247  face_index(const unsigned int i) const;
3248 
3297  neighbor_child_on_subface(const unsigned int face_no,
3298  const unsigned int subface_no) const;
3299 
3325  neighbor(const unsigned int face_no) const;
3326 
3334  int
3335  neighbor_index(const unsigned int face_no) const;
3336 
3344  int
3345  neighbor_level(const unsigned int face_no) const;
3346 
3358  unsigned int
3359  neighbor_of_neighbor(const unsigned int face_no) const;
3360 
3371  bool
3372  neighbor_is_coarser(const unsigned int face_no) const;
3373 
3388  std::pair<unsigned int, unsigned int>
3389  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3390 
3397  unsigned int
3398  neighbor_face_no(const unsigned int neighbor) const;
3399 
3403  static bool
3405 
3419  bool
3420  has_periodic_neighbor(const unsigned int i) const;
3421 
3439  periodic_neighbor(const unsigned int i) const;
3440 
3449  neighbor_or_periodic_neighbor(const unsigned int i) const;
3450 
3466  periodic_neighbor_child_on_subface(const unsigned int face_no,
3467  const unsigned int subface_no) const;
3468 
3479  std::pair<unsigned int, unsigned int>
3480  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3481 
3487  int
3488  periodic_neighbor_index(const unsigned int i) const;
3489 
3495  int
3496  periodic_neighbor_level(const unsigned int i) const;
3497 
3512  unsigned int
3513  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3514 
3520  unsigned int
3521  periodic_neighbor_face_no(const unsigned int i) const;
3522 
3529  bool
3530  periodic_neighbor_is_coarser(const unsigned int i) const;
3531 
3548  bool
3549  at_boundary(const unsigned int i) const;
3550 
3559  bool
3560  at_boundary() const;
3561 
3569  bool
3570  has_boundary_lines() const;
3598 
3616  void
3619 
3623  void
3625 
3633  bool
3635  const unsigned int face_no,
3636  const RefinementCase<dim - 1> &face_refinement_case =
3638 
3644  bool
3645  flag_for_line_refinement(const unsigned int line_no) const;
3646 
3656  subface_case(const unsigned int face_no) const;
3657 
3661  bool
3663 
3668  void
3670 
3674  void
3699  material_id() const;
3700 
3712  void
3713  set_material_id(const types::material_id new_material_id) const;
3714 
3723  void
3724  recursively_set_material_id(const types::material_id new_material_id) const;
3751  subdomain_id() const;
3752 
3768  void
3769  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3770 
3781 
3786  void
3788  const types::subdomain_id new_level_subdomain_id) const;
3789 
3790 
3806  void
3808  const types::subdomain_id new_subdomain_id) const;
3832 
3842 
3856  bool
3857  direction_flag() const;
3858 
3884  unsigned int
3886 
3894  int
3895  parent_index() const;
3896 
3903  parent() const;
3904 
3924  bool
3925  is_active() const;
3926 
3946  bool
3948 
3953  bool
3955 
3989  bool
3990  is_ghost() const;
3991 
3997  bool
3999 
4026  bool
4027  is_artificial() const;
4028 
4035  bool
4037 
4051  bool
4053 
4062  void
4063  set_neighbor(const unsigned int i,
4064  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
4065 
4079  CellId
4080  id() const;
4081 
4083 
4087  double
4088  diameter(const Mapping<dim, spacedim> &mapping) const;
4089 
4107 
4108 protected:
4124  unsigned int
4125  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
4126 
4132  template <int dim_, int spacedim_>
4133  bool
4134  point_inside_codim(const Point<spacedim_> &p) const;
4135 
4136 
4137 
4138 private:
4143  void
4144  set_active_cell_index(const unsigned int active_cell_index) const;
4145 
4149  void
4151 
4155  void
4157 
4161  void
4162  set_parent(const unsigned int parent_index);
4163 
4170  void
4171  set_direction_flag(const bool new_direction_flag) const;
4172 
4173  friend class Triangulation<dim, spacedim>;
4174 
4175  friend class parallel::TriangulationBase<dim, spacedim>;
4176 
4177  friend struct ::internal::TriangulationImplementation::Implementation;
4178  friend struct ::internal::TriangulationImplementation::
4179  ImplementationMixedMesh;
4180 };
4181 
4182 
4183 
4184 /* ----- declaration of explicit specializations and general templates ----- */
4185 
4186 
4187 template <int structdim, int dim, int spacedim>
4188 template <typename OtherAccessor>
4190  const OtherAccessor &)
4191 {
4192  Assert(false,
4193  ExcMessage("You are attempting an illegal conversion between "
4194  "iterator/accessor types. The constructor you call "
4195  "only exists to make certain template constructs "
4196  "easier to write as dimension independent code but "
4197  "the conversion is not valid in the current context."));
4198 }
4199 
4200 
4201 
4202 template <int structdim, int dim, int spacedim>
4203 template <int structdim2, int dim2, int spacedim2>
4206 {
4207  Assert(false,
4208  ExcMessage("You are attempting an illegal conversion between "
4209  "iterator/accessor types. The constructor you call "
4210  "only exists to make certain template constructs "
4211  "easier to write as dimension independent code but "
4212  "the conversion is not valid in the current context."));
4213 }
4214 
4215 
4216 
4217 template <int dim, int spacedim>
4218 template <int structdim2, int dim2, int spacedim2>
4221 {
4222  Assert(false,
4223  ExcMessage("You are attempting an illegal conversion between "
4224  "iterator/accessor types. The constructor you call "
4225  "only exists to make certain template constructs "
4226  "easier to write as dimension independent code but "
4227  "the conversion is not valid in the current context."));
4228 }
4229 
4230 
4231 
4232 template <int structdim, int dim, int spacedim>
4233 template <int structdim2, int dim2, int spacedim2>
4236 {
4237  Assert(false,
4238  ExcMessage("You are attempting an illegal conversion between "
4239  "iterator/accessor types. The constructor you call "
4240  "only exists to make certain template constructs "
4241  "easier to write as dimension independent code but "
4242  "the conversion is not valid in the current context."));
4243 }
4244 
4245 
4246 
4247 template <int dim, int spacedim>
4248 template <int structdim2, int dim2, int spacedim2>
4251 {
4252  Assert(false,
4253  ExcMessage("You are attempting an illegal conversion between "
4254  "iterator/accessor types. The constructor you call "
4255  "only exists to make certain template constructs "
4256  "easier to write as dimension independent code but "
4257  "the conversion is not valid in the current context."));
4258 }
4259 
4260 
4261 #ifndef DOXYGEN
4262 
4263 template <>
4264 bool
4266 template <>
4267 bool
4269 template <>
4270 bool
4272 template <>
4273 bool
4275 template <>
4276 bool
4278 template <>
4279 bool
4281 // -------------------------------------------------------------------
4282 
4283 template <>
4284 void
4286 
4287 #endif // DOXYGEN
4288 
4290 
4291 // include more templates in debug and optimized mode
4292 #include <deal.II/grid/tria_accessor.templates.h>
4293 
4294 #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: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 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
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
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
void set_combined_face_orientation(const unsigned int face, const unsigned char combined_orientation) const
#define DEAL_II_DEPRECATED
Definition: config.h:175
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:475
#define DEAL_II_CXX20_REQUIRES(condition)
Definition: config.h:163
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:476
unsigned int level
Definition: grid_out.cc:4618
#define DeclException0(Exception0)
Definition: exceptions.h:465
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:1614
static ::ExceptionBase & ExcCellFlaggedForRefinement()
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:488
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:510
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: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